Improving the Frequency Precision of Oscillators by Synchronization
IImproving the Frequency Precision of Oscillators bySynchronization
M. C. Cross
Department of Physics, California Institute of Technology, Pasadena CA 91125 (Dated: April 8, 2019)
Abstract
Improving the frequency precision by synchronizing a lattice of oscillators is studied in the phasereduction limit. For the most commonly studied case of purely dissipative phase coupling (theKuramoto model) I confirm that the frequency precision of N entrained oscillators perturbed byindependent noise sources is improved by a factor N as expected from simple averaging arguments.In the presence of reactive coupling, such as will typically be the case for non-dissipatively coupledoscillators based on high-Q resonators, the synchronized state consists of target like waves radiatingfrom a local source which is a region of higher frequency oscillators. In this state all the oscillatorsevolve with the same frequency, however I show that the improvement of the frequency precision isindependent of N for large N , but instead depends on the disorder and reflects the dependence ofthe frequency of the synchronized state on just those oscillators in the source region of the waves. PACS numbers: 05.45.Xt, 87.19.lm, 89.75.Kd a r X i v : . [ n li n . PS ] S e p . INTRODUCTION Oscillators – devices producing a periodic signal at a frequency determined by the char-acteristics of the device and not by an external clock – play a crucial role in much of moderntechnology, for example in timekeeping (quartz crystal watches), communication (frequencyreferences for mixing down radio frequency signals), and sensors. A key characteristic isthe intrinsic frequency precision of the device, which can be quantified in terms of the linewidth of the oscillating signal, or in more detail by the spectral density of the signal in thefrequency domain, or by the Allan deviation in the time domain. This is the fundamentalissue, broadly common to all oscillators, considered in the present paper. There are otherimportant practical characteristics, including the robustness of the frequency to environmen-tal perturbations such as vibrations and temperature fluctuations, that are more dependenton the details of the device implementation; these are not considered here.Unlike a resonator driven by an external oscillating signal, where the line width of thespectral response is determined by the dissipation in the resonator, the line width of anoscillator is only nonzero in the presence of noise. An oscillator is mathematically describedby a limit cycle in the phase space of dynamical variables, and the line width of the signalcorresponding to a limit cycle is zero. Dissipation serves to relax the system to the limitcycle, which is itself determined by a balance of energy injection and dissipation, and doesnot broaden the spectral line. The spectral line is broadened only if there is some stochasticinfluence that causes the phase space trajectory to fluctuate away from the limit cycle. Thusthe frequency imprecision (line width) is due to noise[1].One way to improve the frequency precision of oscillators that has been suggested invarious scientific disciplines[2, 3], is to sum the signal from a number N similar oscillators: inthis case, if the noise sources are uncorrelated over the individual oscillators, simple averagingsuggests that the effective noise intensity will be reduced by a factor 1 /N . Of course, dueto fabrication imperfection the isolated oscillators will have slightly different frequencies,and so simply averaging the summed signal from the individual oscillators will also tendto average the signal to zero. In addition, the line width of the reduced intensity signalwill reflect the frequency dispersion of the devices. If however some coupling is introducedbetween the oscillators, they may become synchronized to a state in which all the oscillatorsare entrained to run at a single common frequency[4]. In this case, the averaging argument2ould suggest a 1 /N reduction in the effective noise, and so a factor of N enhancement inthe frequency precision.In this paper I investigate the improvement of frequency precision due to synchronizationin canonical models of the phenomenon. I focus on lattices of oscillators with nearest neigh-bor coupling. I confirm that the factor N improvement applies exactly in the case of purelydissipative coupling between the phases of the oscillators. However, this scaling breaks downif there is in addition a reactive (non-dissipative) component of the coupling. Examples ofreactive coupling are a displacement, rather than velocity, coupling for arrays of mechanicaloscillators[5, 6] or trapped ions[7]. In this case the improvement factor becomes independentof N for large N , but depends on the amount of disorder. This may lead to much poorerfrequency precision than anticipated from the na¨ıve averaging argument.The main focus of the paper is on d-dimensional lattices of oscillators with nearest neigh-bor coupling. I briefly discuss the extension to longer range coupling and to complexnetworks. These results should be relevant to questions of the precision of synchronizedoscillations in biological contexts. II. MODEL
The model I consider is N nearest neighbor coupled phase oscillators[4, 8–11]˙ θ i = ω i + (cid:88) nn j Γ( θ j − θ i ) + ξ i ( t ) , i = 1 , . . . , N. (1)Here ω i are the frequencies of the individual oscillators, which are assumed to be independentrandom variables taken from a distribution g ( ω ) with width σ [12]. The term ξ i ( t ) representsthe noise acting on the i th oscillator. I will assume white noise, but the results are easilygeneralized to colored noise, and also to noise that depends on the phase θ i . An importantassumption is that the noise is uncorrelated between different oscillators. Thus (cid:104) ξ i ( t ) ξ j ( t (cid:48) ) (cid:105) = c ( t − t (cid:48) ) δ ij , (2)where for white noise c ( t − t (cid:48) ) = f δ ( t − t (cid:48) ) , (3)with f the individual oscillator noise strength, taken to be the same for all oscillators, sinceI am imagining a system where the oscillators are designed to be as similar as possible.3he coupling between the oscillators is given by the function Γ, a 2 π periodic function ofthe phase differences of nearest neighbor oscillators. A commonly used model[8, 9] is givenby the coupling function Γ( φ ) = sin φ. (4)Any parameter K multiplying sin φ and giving the strength of the coupling may be scaledto unity by rescaling time and frequencies. Thus the only parameters defining the behaviorof the system are the distribution of the frequencies ω i , and in particular the width σ of thedistribution after this rescaling (i.e., the width of the frequency distribution relative to thecoupling strength).The coupling function Eq. (4) is antisymmetric Γ( − φ ) = − Γ( φ ), and equation (1) ispurely dissipative[13]. A more general coupling functionΓ( φ ) = sin φ + γ (1 − cos φ ) (5)breaks this symmetry for nonzero γ , and includes non-dissipative, propagating effects[10, 14].I will use this model to study the effect of reactive coupling. Without loss of generality, Itake γ > φ ) (cid:39) γ − ( e γφ − , | φ | , | γφ | (cid:28) . (6) III. SYNCHRONIZATION
I first describe the behavior predicted by Eqs. (1-5) in the absence of noise. For sufficientlyweak disorder (small σ ) and for a finite number of oscillators, the oscillations described byEqs. (1-3) with f = 0 become entrained in the sense that all the phases advance at the sameconstant rate ˙ θ i = Ω . (7)The solution to these equations can be written θ i ( t ) = θ ( s ) i + Θ( t ) , (8)4ith θ ( s ) i a fixed point solution in the rotating frame, and Θ the phase of the collective limitcycle given by Θ( t ) = Ω t + Θ , (9)with Θ an arbitrary constant. θ ( s ) i and Ω are given by solvingΩ = ω i + (cid:88) nn j Γ( θ ( s ) j − θ ( s ) i ) , (10)when solutions exist. The behavior in the thermodynamic limit N → ∞ for differentlattice dimensions d , and the critical value of σ for the onset of this entrained state andits dependence on system size, lattice dimension, frequency distribution etc., are subtlequestions that have not been fully answered. But for the practical case of a finite number ofoscillators, we expect that an entrained (fully frequency locked) state will exist for sufficientlysmall σ [16]. Such a state is a limit cycle of the system of oscillators with frequency Ω. Bymoving to a rotating frame, θ i ( t ) → θ i ( t ) − Ω t , the entrained state becomes a fixed point,simplifying the subsequent analysis.The nature of the synchronized state depends sensitively on whether the coupling is purelydissipative ( γ = 0) or also contains a reactive component ( γ (cid:54) = 0). For purely dissipativecoupling the interactions cancel when summed over a block of oscillators. This means thatthe frequency of the entrained state is the mean ¯ ω of the oscillator frequencies. For a one-dimensional lattice, the individual phases θ ( s ) i are then given by θ ( s ) i +1 − θ ( s ) i = − sin − X i with X i = (cid:80) ij =1 ( ω j − ¯ ω )[17]. The accumulated randomness X i performs a random walk asa function of the lattice index i , and so the strain θ ( s ) i +1 − θ ( s ) i also varies with i roughly asa random walk. The break down of the synchronized state as the disorder or system sizeincreases occurs when the excursion of X i exceeds unity (remember the coupling strengthis scaled to one). This occurs for σ ∼ N − / . On the other hand, for γ (cid:54) = 0, the interactionterms summed over a block of oscillators do not cancel, and the same arguments cannot bemade. For γ > λ propagating away from a unique source in the system, located ata cluster of higher frequency oscillators[14, 15]. The derivation of this result is described inmore detail below. In this wave state the frequencies of all the oscillators remain entrained,even though the phases vary by more than 2 π over the system for λ < N , and by manyfactors of 2 π for λ (cid:28) N . In two dimensional lattices, roughly circular “target” waves are5ound for γ (cid:54) = 0, with the waves propagating away from a source with location again givenby a core region of higher frequency oscillators[14, 15]. IV. FREQUENCY PRECISION AND NOISE
The noise terms in Eqs. (1) will lead to deviations of the solution from the limit cycleEqs. (8,9) and so to a broadening of the spectral lines of the output signal from the oscillator.The full noise spectrum depends on a complete solution of Eqs. (1). However, for frequencyoffsets from the no-noise peaks in the power spectrum that are small compared with therelaxation rates onto the limit cycle, the effects of the noise can be reduced to a singlestochastic equation for the limit cycle phase Θ that gives the collective behavior of theentrained oscillators. This result has been derived for a general limit cycle by a numberof authors using a variety of formalisms[1, 18–20]. The key idea is that a change in Θcorresponds to a time translation, and so gives an equally good limit cycle solution: thus aperturbation to Θ does not decay, and this represents a zero-eigenvalue mode of the linearstability analysis[21]. The remaining eigenvalues of the stability analysis will be negative,corresponding to exponential decay onto the limit cycle. For time scales longer than theserelaxation times, it is only the projection of the noise along the zero eigenvalue eigenvector,that is important: the other fluctuation components will have decayed away.For the white noise sources considered here, this stochastic equation for the phase issimply[20] ˙Θ( t ) = ¯Ω + Ξ( t ) , (11)with (cid:104) Ξ( t )Ξ( t (cid:48) ) (cid:105) = F δ ( t − t (cid:48) ) , (12)with F the noise strength resulting from the projection and ¯Ω the limit cycle frequencywhich is Ω with an O ( F ) correction. The solution to Eq. (11) is a drift of the mean phaseat the rate ¯Ω (cid:104) Θ( t ) (cid:105) = ¯Ω t, (13)together with phase diffusion (cid:104) (Θ( t ) − ¯Ω t ) (cid:105) = F t. (14)6n output signal from the oscillator such as X = cos Θ( t ) will have a power spectrumconsisting of a Lorentzian peak centered at ¯Ω (and, for more general signals, the harmonics) S XX ( ω ) = S π F ( ω − ¯Ω) + F , (15)with S the spectral weight of the delta-function peak in the spectrum of the no-noiseoscillator[19]. The width of the spectral peak is therefore equal to the phase noise strength F . Thus the tails of the spectrum away from the peaks decay as ω − ; this is the white-noisecomponent of the Leeson noise spectrum for oscillators[22]. Other noise spectra will lead todifferent power law tails.The relationship of the effective noise strength F acting on the collective phase Θ tothe the strength of the noise f acting on each individual oscillator is given by projectingthe individual noise components ξ i ( t ) along the phase variable Θ. Denoting the phases θ i by the vector θ , the tangent vector to the limit cycle (the zero-eigenvalue eigenvector) isgiven by e = (1 , , . . . , t are δθ i = e ,i ∆ t . Using the general results of refs. 18–20, or the simpler analysis for the present case sketched in Appendix A, the relationshipis F = e † · e † ( e † · e ) f, (16)with e † the zero-eigenvalue adjoint eigenvector. Thus finding the broadening of the line dueto the noise is reduced to calculating the adjoint eigenvector e † .The Jacobean matrix J yielding the linear stability analysis of the phase dynamics aboutthe fixed point phases θ ( s ) i defining the limit cycle is J ij = Γ (cid:48) ( θ ( s ) j − θ ( s ) i ) for i,j nearest neighbors , (17a) J ii = − (cid:88) nn j Γ (cid:48) ( θ ( s ) j − θ ( s ) i ) , (17b)with other elements zero. The vectors e , e † are defined by J · e = 0 , (17c) J † · e † = 0 , (17d)with J † ij = J ji . 7ote that I am treating noise perturbatively in the small noise limit and for a finite system:in this case the result is given by just the effect on the overall phase of the synchronizedstate, which is the zero-mode of the system. I am not considering modifications to thesynchronized state due to the noise, such as changes in values of the critical disorder forsynchronization, or changes in the nature of the synchronized state. In a finite system therewill be barriers to such fluctuations, and their rates will vary with the noise strength f as e − ∆ /f with ∆ some number depending on the states considered. These fluctuations cantherefore be ignored for small enough f . As the number of oscillators tends to infinity, somebarriers will become very small, and the synchronized state may be significantly changed oreven eliminated by the addition of noise[11], as for phase transitions in equilibrium systemsat finite temperature. V. DISSIPATIVE COUPLING
For purely dissipative coupling γ = 0, the Jacobean J is symmetric and so the adjointeigenvector is equal to the forward eigenvector which is the tangent vector defined by aninfinitesimal time translation e † = e = (1 , , . . . , . (18)This result is true for any antisymmetric coupling, and is not restricted to the nearestneighbor model. This immediately gives the result for the effective noise strength F = N − f, (19)so that the frequency precision of the entrained state is enhanced by the factor N . Notethat the enhancement does not depend on the degree of phase alignment quantified by themagnitude of the order parameter Ψ = N − (cid:80) j e iθ ( s ) j , which may be less than unity (i.e.,phases not fully aligned) even in the entrained state. The result Eq. (19) has been obtainedpreviously[2, 3], although I believe the present derivation is more systematic, since it doesnot assume that the effect of the noise on the collective phase remains small at long times.8 I. DISSIPATIVE PLUS REACTIVE COUPLING
For general coupling the Jacobean is not symmetric, and there is no obvious relationshipbetween e † and e in general. Physically, in situations where the entrained state consistsof waves emanating from a source region of higher frequency oscillators, we might expectthe frequency precision to be determined by fluctuations of only those oscillators in the coreregion that fix the frequency of the waves. This means that the reduction of the effectivenoise by averaging is only over this core region of oscillators, giving a poorer improvementof the frequency precision. I first demonstrate this result for a simpler “one-way” couplingfunction introduced by Blasius and Tonjes[15] for which analytic solution is possible. I thenderive the result for the general coupling function Eq. (5) assuming the disorder is smallenough so that the phase difference between all nearest neighbor oscillators is small, in whichcase the approximation Eq. (6) may be used. The result depends on the mapping[14, 15]of the solution for the entrained state onto the Anderson localization problem[23], and theknown properties of the localized states in this problem[24, 25], together with a relationshipbetween e † and the localized states that I demonstrate. I also investigate one and twodimensional lattices numerically. A. One-way coupling
Blasius and Tonjes[15] proposed a simple, exactly soluble model of a one dimensionallattice, with a nearest neighbor coupling function such that the phase of oscillator i is onlyinfluenced by its neighbors if their phases are ahead (all phase differences are assumed to besmall so that this notion makes sense). I use the exampleΓ( φ ) = γ ( e γφ −
1) for φ > , φ < . (20)The entrained solution is given by Ω = ω m with m the index of the largest frequency in thelattice, and then the fixed point solution θ ( s ) is constructed iteratively from θ ( s ) m using θ ( s ) i = θ ( s ) i − − γ − ln[1 + γ ( ω m − ω i )] for i > m,θ ( s ) i +1 − γ − ln[1 + γ ( ω m − ω i )] for i < m. (21)The value chosen for θ ( s ) m sets the overall phase Θ. An example of the entrained state for200 oscillators in a 1d lattice is given in Fig. 1, showing waves emanating from the oscillator9 (cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230) (cid:45) (cid:45) (cid:45) (cid:45) (cid:45)
50 Index P h a s e Θ FIG. 1. Entrained state for the one sided model Eq. (20) for 200 oscillators. Left panel: oscillatorphase as a function of lattice site with θ m set to 0; right panel: gray scale plot of cos θ i as a functionof time. The oscillator frequencies were taken from a uniform frequency distribution with width σ = 0 . γ = 1. with maximum frequency at m = 121.It is easy to see for this coupling function that the zero-eigenvalue adjoint eigenvector is e † ,i = δ im , (22)corresponding to the fact the phase θ m is not coupled to either neighbor, since θ m > θ m ± .Thus the effective noise is given by F = f , and there is no improvement of the frequencyprecision, even though all the oscillators are entrained. This is because the single oscillatorwith maximum frequency determines the entrained frequency, and therefore the entrainedfrequency is as sensitive to noise as this single oscillator. B. General coupling
I now consider the case of general coupling Eq. (5) in the limit of small enough disorderso that the small phase difference approximation Eq. (6) approximation may be used. Inthis case, Blasius and Tonjes[15] showed that the Cole-Hopf transformation θ i = γ − ln q i q i = Eq i = γω i q i + (cid:88) nn j ( q j − q i ) , (23)with the eigenvalue E = γ Ω. This is equivalent to the tight binding model for a quantumparticle on a random lattice, and the properties of the solution can be extracted in analogywith Anderson localization[23]. At long times the solution q ( t ) = q max e E max t correspondingto the largest eigenvalue E max will dominate. This gives the entrained state θ ( s ) i = γ − ln q max i , (24)with frequency Ω = E max /γ . Anderson localization theory shows that q max may be chosenpositive, and it has the form of an exponentially localized state centered on a region of thelattice with a concentration of larger frequency oscillators. The exponential localization of q max corresponds to a roughy linear phase profile, again leading to waves propagating froma source, as shown in Fig. 2.I now analyze the frequency precision based on the properties of the solution q max knownfrom studies of the Anderson problem.
1. One-dimensional lattice
For a one-dimensional lattice with nearest neighbor coupling Eq. (6), the Jacobean matrixEq.˜(17) for the stability analysis of the fixed point solution θ ( s ) i is the tridiagonal matrix withelements J ii ± = e γ ( θ ( s ) i ± − θ ( s ) i ) = q max i ± /q max i , (25) J ii = − J i,i +1 − J i,i − , (26)except for the first and last rows corresponding to the end oscillators which only have oneneighbor J = − J = q max2 /q max1 , J NN − = − J NN = q max N − /q max N , (27)and all other elements zero. It is easily checked that (1 , , . . . ,
1) is indeed the zero-eigenvalue eigenvector. The adjoint matrix has off diagonal elements J † ii ± = q max i /q max i ± (28)11 (cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230) (cid:45) (cid:45) (cid:45) (cid:45) (cid:45) (cid:45) P h a s e Θ FIG. 2. Entrained state for the general model Eq. (5) for a chain of 200 oscillators with nearestneighbor coupling, using the small phase difference approximation Eq. (6). Left panel: oscillatorphase as a function of lattice site; right panel: gray scale plot of cos θ i as a function of time. Theoscillator frequencies were taken from a uniform frequency distribution with width σ = 0 . γ = 1. except for the first and last rows for which J † = q max1 /q max2 , J † NN − = q max N /q max N − , (29)and diagonal elements J † ii = J ii , (30)with all other elements zero. The key result is that the (unnormalized) adjoint eigenvectorcan be found explicitly e † ,i = ( q max i ) , (31)as can be confirmed by direct substitution. This simple result follows from the quotientform of the Jacobean matrix elements for the special form of the interaction Eq. (6). Anexample of the vector q max and the adjoint eigenvector e † for the system of Fig. 2 is shownin Fig. 3.The noise reduction factor F/f , Eq. (16), is given by Ff = (cid:80) i ( q max i ) (cid:2)(cid:80) i ( q max i ) (cid:3) . (32)12 (cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:230)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)(cid:224)
80 100 120 1400.00.20.40.60.81.0 Index q , e γσ N o i s e r edu c t i on f a c t o r FIG. 3. Left panel: Localized solution q max (blue circles) and zero-eigenvalue adjoint eigenvector e † (black squares) for the system of Fig. 2. (The normalizations are chosen for the plot so thatthe largest element of each vector is 1.) Only the portion of the L = 200 lattice where the elementshave appreciable size is shown. Right panel: scaling of the noise reduction factor F/f deducedfrom Eq. (32) with γσ . Each point is the average of 1000 realizations of the random lattice offrequencies. This equation directly relates the improvement in frequency precision to the solution ofthe linear Anderson problem Eq. (23). The expression Eq. (32) for the noise reduction isthe inverse participation ratio p − of the vector q max of the linear localization problem.This can be used to define the radius of the localized state r ≡ p/
2. Thus I find that thenoise reduction factor is given by the size of the source of the waves, rather than by thetotal number of oscillators, giving an improvement in frequency stability that is significantlyworse for a large number of entrained oscillators. The size of the source is defined precisely interms of the participation ratio of the maximum energy localized state of the correspondingAnderson problem. For the system in Figs. 2,3, p (cid:39) .
36. In this example, the frequencyprecision would not be improved by increasing the number of oscillators beyond about ten.From Eq. (23) it is clear that the noise reduction factor
F/f depends on the parametersof the model only through the product γσ , for the approximation to Γ( φ ) used. Within thisapproximation, the scaling found from numerical solutions of Eq. (23) for a one dimensional13attice is shown in Fig. 3. The calculations were done for chains of length 100-1000, andthe results were insensitive to the length providing it is much larger than the width of thelocalized state. A power law F/f ∝ ( γσ ) . is a good fit to the calculated results over therange considered.
2. d-dimensional lattice
FIG. 4. Source and waves for a 60 x 60 two dimensional lattice of oscillators with γ = 1 , σ = 2.Left panel: zero-eigenvalue adjoint eigenvector e † showing the effective size of the source of thewaves in the entrained state; right panel: grey scale plot of cos θ ( s ) i,j giving a snapshot of the wavesemanating from the source. The same argument applies to general dimension, although the structure of the Jacobeanmatrix is no longer tridiagonal. Choose any convenient labeling of the oscillators θ i , i =1 . . . N . Nearest neighbor oscillators will not in general be adjacent in the list. However theJacobean and its adjoint are still defined by J ij = q max j /q max i , for ij nearest neighbors (33) J † ij = q max i /q max j , for ij nearest neighbors (34) J † ii = J ii = − (cid:88) nn j q max j /q max i , (35)with other elements zero. The eigenvectors e , e † are as before, and the expression Eq. (32)for F in terms of the inverse participation ratio is unchanged. Thus I expect the noise14eduction factor to scale as r − d with r ∼ p /d the radius of the maximum energy localizedstate in the d -dimensional Anderson localization problem.Figure 4 shows an example of a target wave entrained state for a 60 x 60 two dimensionallattice. The left panel is the adjoint eigenvector defining the source: the inverse participationratio, yielding the improvement in the frequency precision, is 54.1 (cf. N = 3600). The rightpanel is a plot of cos θ ( s ) ij , calculated from the Cole-Hopf transformation and the eigenvector q max , giving a snapshot of the waves in the entrained state. C. More general systems
The results Eq. (16) and Eq. (31) remain valid for a more general coupling yielding theequations for the phase dynamics˙ θ i = ω i + (cid:88) j K ij Γ( θ j − θ i ) + ξ i ( t ) , i = 1 , . . . , N. (36)with K ij a symmetric matrix giving the strength of the coupling between oscillators i and j . The small phase difference condition so that Eq. (6) may be used is now that the phasedifference between any two oscillators with nonzero K ij be small in the entrained state. Thesame analysis leads to the relationship Eq. (32) between the noise reduction factor and theparticipation ratio of the largest E eigenvector q max of the corresponding linear problem Eq i = γω i q i + (cid:88) j K ij ( q j − q i ) . (37)Note that although the strength of the coupling can be different for different paris of oscil-lators, the form of the coupling Eq. (5), and in particular the ratio of reactive to dissipativecomponents, must be the same for this simple analysis to apply.One generalization Eq. (36) allows is to a lattice of oscillators with short range, butnot just nearest neighbor, interactions. The scaling of the noise reduction with γσ will bethe same as for nearest neighbor interactions, since the scaling properties of the Andersonproblem are the same for these two cases. More generally, the method reduces the problemof calculating the improved frequency precision in the entrained state of a complex networkof oscillators to solving the linear problem Eq. (37) for the network architecture and couplingparameters K ij . 15 II. DISCUSSION
The major result of this paper is that for oscillators on a lattice with short range couplingincluding a reactive component, the improvement of the frequency precision due to synchro-nization is limited to a factor given by the number of oscillators in the core source region ofthe waves that form the entrained state, rather than a factor equal to the total number ofoscillators, as is the case for purely dissipative coupling. I showed this result explicitly forthe phase reduction description, in the limit of small enough disorder, or strong enough cou-pling, so that the phase differences between interacting oscillators are small in the entrainedstate. The size of the core region is given by the extent of the localized ground state of thecorresponding linear Anderson problem, onto which the nonlinear phase equation is mappedby a Cole-Hopf transformation. The precise relationship is Eq. (32) relating the reductionin the phase noise to the inverse participation ratio of the localized state. This relationshipremains true for general networks of oscillators providing the small phase difference approx-imation Eq. (6) applies for all interacting pairs of oscillators, and reduces the calculation ofthe frequency precision to the corresponding linear Anderson problem on the network.Within the small phase difference approximation, the entrained state of waves propagatingfrom the localized source is the unique state at long times. However, for the phase equationswith the full coupling function Eq. (5) other states may result depending on the initialconditions. This is particularly evident for two dimensional lattices, where spiral statesare seen in numerical simulations starting from particular initial conditions[14]. Due tothe topological constraint of integral 2 π phase winding around the center, such a structuresurvives at long times, unless the core migrates to an open boundary. A second consequenceof the topological structure is that there are necessarily large phase differences betweennearest neighbor phases in the core, so that the small phase difference approximation breaksdown. In the spiral state all the oscillators are again entrained to a single frequency, thatprobably depends just on the oscillators in the core region. It would be interesting to extendthe analysis of the present paper to these spiral states.16 ppendix A: Derivation of the phase equation In this appendix I give a brief derivation of the stochastic phase equation Eq. (11). Thisequation follows from the general results for limit cycles of refs. 19 and 20, but the derivationis simpler for the phase reduction description. I start from Eq. (1)˙ θ i = ω i + (cid:88) nn j Γ( θ j − θ i ) + ηξ i ( t ) , i = 1 , . . . , N, (A1)introducing a perturbation parameter η to label the small noise. I will develop a perturbationexpansion in η , and set η → η to extract the phase diffusion: continuing to second order would be needed to find theLamb shift like correction to the mean frequency.At zeroth order in η the solution is the no-noise solution Eqs. (8-10), with Θ an arbitraryconstant. In the presence of order η noise I expect this overall phase to evolve on a slowtime scale T = ηt , Θ = Θ ( T ), so that the phases will be given up to order η by θ i ( t ) = θ ( s ) i + Ω t + Θ ( T ) + ηθ (1) i ( t ) + . . . , (A2)with θ (1) ( t ) a correction to be found. Expanding the equation of motion up to order η theequation for this correction is˙ θ (1) i − (cid:88) j J ij θ (1) j = − (Θ (cid:48) e ,i − ξ i ( t )) , (A3)with J ij the Jacobean Eq. (17), Θ (cid:48) = d Θ /dT , and e = (1 , , . . . θ (1) along eigenvectors of J with negative eigenvalues will have some finite value given byinverting this equation. However, for the component along the zero-eigenvalue eigenvector,there is no restoring force, and any nonzero value of the right hand side will lead to largevalues of θ (1) i at large times, violating the assumption that θ (1) gives a small correction.This component is extracted by multiplying on the left by the adjoint eigenvector e † since e † · J = 0. This leads to the solvability condition for θ (1) to remain finiteΘ (cid:48) = e † · ξ e † · e . (A4)Returning to the original variables and setting η → t ) = Ω + Ξ( t ) , (A5)17ith the correlation function of the effective noise (cid:104) Ξ( t )Ξ( t (cid:48) ) (cid:105) = (cid:80) ij e † ,i e † ,j (cid:104) ξ i ( t ) ξ j ( t (cid:48) ) (cid:105) ( e † · e ) , (A6)giving Eq. (12) with Eq. (16) for equal, uncorrelated white noise of strength f acting on eachindividual oscillator. Note that the result does not depend on the choice of normalizationfor e † . A specific normalization choice for e was made in setting up Eq. (A3). ACKNOWLEDGMENTS
This work was supported by NSF grant number DMR-1003337. I thank Eyal Kenig andTony Lee for comments on the manuscript. [1] M. Lax, Phys. Rev. , 290 (1967).[2] H. C. Chang, X. D. Cao, U. K. Mishra, and R. A. York, IEEE Transactions on MicrowaveTheory and Techniques , 604 (1997).[3] D. J. Needleman, P. H. E. Tiesinga, and T. J. Sejnowski, Physica D , 324 (2001).[4] A. Pikovsky, M. Rosenblum, and J. Kurths, Synchronization: A universal concept in nonlinearscience (Cambridge University Press, Cambridge, 2001).[5] M. C. Cross, A. Zumdieck, R. Lifshitz, and J. L. Rogers, Phys. Rev. Lett , 224101 (2004).[6] M. C. Cross, J. L. Rogers, R. Lifshitz, and A. Zumdieck, Phys. Rev. E , 036205 (2006).[7] T. E. Lee and M. C. Cross, Phys. Rev. Lett. , 143001 (2011).[8] A. T. Winfree, J. Theor. Bio. , 15 (1967).[9] Y. Kuramoto, Lecture Notes in Physics , 420 (1975).[10] Y. Kuramoto, Chemical Oscillations, Waves, and Turbulence (Springer, New York, 1984).[11] J. A. Acebron, L. L. Bonilla, C. J. P. Vicente, F. Ritort, and R. Spigler, Rev. Mod. Phys. , 137 (2005).[12] For distributions with finite variance, σ can be defined as the standard deviation. For distri-butions without a variance, such as a Lorentzian, some other appropriate characterization ofthe width may be chosen.[13] D. Topaj and A. Pikovsky, Physica D , 119 (2002).
14] H. Sakaguchi, S. Shinomoto, and Y. Kuramoto, Progress of theoretical physics , 1069(1988).[15] B. Blasius and R. Tonjes, Phys. Rev. Lett. , 084101 (2005).[16] There may be multistability, so that other non-entrained states may exist for the same pa-rameter values: I will only consider the entrained state.[17] S. H. Strogatz and R. E. Mirollo, Physica D , 143 (1988).[18] F. X. Kaertner, Int. J. Circuit Theory Applicat. , 485 (1990).[19] A. M. A. Demir and J. Roychowdhury, IEEE Trans. Circuits and Syst. , 655 (2000).[20] D. S. Goldobin, J.-N. Teramae, H. Nakao, and G. B. Ermentrout, Phys. Rev. Lett. (2010).[21] In the general case the stability analysis would be a Floquet analysis of a periodic state; inthe present case this can be reduced to a stability analysis of the fixed point in the rotatingframe.[22] D. B. Leeson, Proc. IEEE , 329 (1966).[23] P. W. Anderson, Phys. Rev. , 1492 (1958).[24] P. A. Lee and T. V. Ramakrishnan, Rev. Mod. Phys. , 287 (1985).[25] M. Jensen, Phys. Rep. , 2 (1998)., 2 (1998).