Characteristic Lyapunov vectors in chaotic time-delayed systems

0.53 MB·

by  Diego Pazó

submittedtoPhys.Rev.E
Characteristic Lyapunov vectors inchaotic time-delayedsystems
Diego Pazo´ and Juan M. Lo´pez
InstitutodeF´ısicadeCantabria(IFCA),CSIC–Universidad deCantabria, E-39005Santander, Spain
(Dated:January17,2011)
WecomputeLyapunovvectors(LVs)correspondingtothelargestLyapunovexponentsindelay-differential
equationswithlargetimedelay. WefindthatcharacteristicLVs,andbackward(Gram-Schmidt)LVs,exhibit
long-rangecorrelations,identicaltothosealreadyobservedindissipativeextendedsystems.Inadditionwegive
numericalandtheoreticalsupporttothehypothesisthatthemainLVbelongs,underasuitabletransformation,
totheuniversalityclassoftheKardar-Parisi-Zhangequation. Thesefactsindicatethatinthelargedelaylimit
(animportantclassof)delayedequationsbehaveexactlyasdissipativesystemswithspatiotemporalchaos.
PACSnumbers:05.45.Jn,02.30.Ks,05.40.-a
1
1
0
I. INTRODUCTION delay:
2

n dy
a Delayeddynamicalsystems(DDSs)servetomodeldiverse =F(y,yτ) (1)
dt
J phenomena in physics [1] (prominently in optics), but also
4 in other fields [2] such as engineering, biology, climatology wherey =y(t−τ). InaLyapunovanalysisweareinterested
τ
1 or ecology. More than one decade ago [3–5] it was found intheevolutionofinfinitesimalperturbationsδy(t):
that delayed systems with one constantdelay can be studied

D] like extended systems, and they present not only analogies dδy
butequivalentphenomenasuchaspatterninstabilities[6]and =uδy+vδyτ, (2)
dt
C
high-dimensionalchaos[7]. ChaosinDDSsisconsideredas
n. aparticulartypeofspatio-temporalchaosforwhichthedelay whichisalso a delayedequationwithδyτ ≡ δy(t−τ), and
i playstheroleofthesystemsize. where u and v are functions: u(y,yτ) ≡ ∂yF, v(y,yτ) ≡
l
n There are a number of studies concerning the tangent dy- ∂yτF. If the DDS is chaotic an initial random perturbation
[ namics of systems with time delay. Previous works mainly becomesaftersometimealignedwiththemainLV(reaching

focusedontheLyapunovexponents(LEs)andpropertiesob- astationarystateinastatisticalsense),andtheaverageexpo-
1
tained fromthem(dimension,entropy,...) [5, 7–10]. Less is nentialgrowthrateisthelargestLEofthesystemλ . Numer-
v 1
9 known aboutthe associated tangent space directions, in par- ically,thisissimplyachievedbyintegratingEqs.(1) and(2)
7 ticular, there are no studies analyzing the structure of Lya- for long enoughtimes. For non-leadingLVs, more involved
7 punovvectors(LVs)apartfromthemainone(pointingalong algorithmsareneeded(seebelow).
2 the most unstable direction). Correlations of the main LV Wefocusourstudyonsystemsinwhichdelayedandnon-
.
1 wereonlyrecentlyaddressed[11,12]takingadvantageofthe delayedtermsareseparated:F(y,y )=Q(y)+R(y ). Our
τ τ
0 genericmappingofaDDSintoanextendedsystem,butwith resultsareexpectedtobegenericforthisclassofmodels. In
1 contradictingresults. our study, we have assumed the non-delayed part is linear1,
1
When consideringLVs otherthan the main one, one must Q(y) = −ay, and all nonlinearities appear in the retarded
:
v distinguishamongdifferentvectortypes. Theso-calledchar- component R. Many important time-delayed systems, in-
i acteristicLVs [13]constitutetheonlyintrinsicsetofvectors cludingtheMackey-Glass(MG)[15]andIkedamodels[16],
X
that is univocally defined and is covariant with the dynam- andopticaldelayedfeedbacksystems[17],canbeexpressed
r
a ics. However,theirnumericalcomputationisdifficultinhigh- in this mathematical form. In particular, we have studied
dimensionalsystems. We adapttoDDSsoneofthemethods numericallydifferentnonlinearfunctions: the Mackey-Glass
tocomputecharacteristicLVs[14]. (MG),R(ρ)=bρ/(1+ρ10),[15];andthenonlinearfunction
In this paper we reporton the generic propertiesthat LVs R(ρ)=bsin2(ρ−ψ)thatappearsinsomeopticalcryptosys-
exhibitin DDSs in the large delay limit. Our numericaland temswithdelayedfeedback[17,18]. Almostidenticalresults
theoretical results indicate that the LVs (the main one and inqualitativeandquantitativetermsareobtainedforthesetwo
the others) behave in qualitative and quantitative terms like systems. Therefore, for the sake of brevity, we choose to
in one-dimensional dissipative systems with spatiotemporal presenttheresultsonlyfortheMGmodel,alsousedin[11].
chaos. In our simulationswe have integratednumericallyEq. (1)
using a third-order Adams-Bashforth-Moulton predictor-
corrector method [19], while linear equations, e.g. (2), have
been integratedusing an Euler methodwith the non-delayed
II. TIME-DELAYEDCHAOTICSYSTEMS
DDSs may exhibit chaos even in the simplest situation in
which the main variable is a scalar and its evolution is de- 1WecarriedoutsomesimulationswithQ(y) = −ay3/2 thatdidnotre-
termined by a delay-differential equation with one constant vealedrelevantdifferences.
2
part integrated semi-implicitly. The parameters we used for LV one has to compute the first n−1 forward LVs in addi-
the MG model are a = 0.1, b = 0.2 and the time step was tiontothefirstnbackwardLVs. ForwardLVs{f (x,θ)}are
n
dt = 0.2. In a numericalintegration of a time-delayed sys- obtainedlikebackwardLVsbutintegratingtheperturbations
tem [7], the temporaldiscretization makes the system finite- backwardintime(fromtheremotefuturetothepresenttime
dimensional (with τ/dt degrees of freedom) but this has no t). They have to be calculated using the transpose Jacobian
effectintheresultswhatsoeversinceanincreaseoftimeres- matrixsothatLEsareobtainedwiththeusualordering[21]:
olution (i.e. decreasing dt) does not modify the largest LEs λ ≥ λ ≥ ···. Forourtime-delayedsystemthisprocedure
1 2
(andtheirassociatedvectors). encompassestointegrate:
In several previousworks [3–5, 11, 12] it was found use-
fultomaptheDDSintoanequivalentspatiallyextendeddy- dδy(t)
namicalsystemof“size”τ thatevolvesatdiscretetimesθ as − =uδy(t)+v˜δy(t+τ). (4)
dt
follows.Weexpressthecontinuoustimeas
The diagonalcoefficient u = ∂ F(w,z)| is
t=x+θτ (3) w w=y(t),z=y(t τ)
identicaltothatinEq.(2).Theoff-diagonaltermhasacu−rious
wherex∈[0,τ)isthe“spatial”positionandθ ∈Zisthedis- structure stemming from transposing a Jacobian matrix with
cretetime. Thisspatialrepresentationcanbeappliedtoboth, nonzeroelementsτ temporalunitsbelowthe maindiagonal.
thestateofthesystem,y(t)→y(x,θ),andtheperturbations, Thus, in Eq. (4), we have v˜ = ∂zF(w,z)|w=y(t+τ),z=y(t).
δy(t)→δy(x,θ). Thebenefitofthisspatialrepresentationis Notice that, as the integration runs backward in time (t →
thatonecananalyzetheDDSwiththetoolsavailableforspa- −∞), the minus sign in the left hand side is canceled. As
tially extended dynamical systems. We shall use the spatial occurswithbackwardLVs,periodicGram-Schmidtorthonor-
representationofdelayedsystemsheretoanalyzeLVsandan- malizations are needed to avoid the collapse of all perturba-
alyzetheexistenceoflong-rangecorrelationsaswellassome tionsalongthemostunstabledirection.Notealsothat,dueto
otherdynamicalproperties. thedelay,Eq.(1)cannotbeintegratedbackward,andhenceit
is necessary to store the trajectory y(t) at every time step in
thecomputertobeusedinthecomputationof(4).
III. CHARACTERISTICVS.BACKWARDLYAPUNOV OncebothbackwardLVsandforwardLVshavebeencom-
VECTORS puted the characteristic LVs are easily calculated following
the prescriptionsin Ref. [14]. One delicate point in delayed
Thereissomedegreeofconfusionintheliteratureregard- systems is that one must make sure that both forward and
ing the definition and computation of LVs. Many authors backwardLVsetscorrespondexactlytothesametimeinterval
thinkofLVsastheorthonormalframeofvectorsthatresults [t,t+τ).
as a byproduct of computing the LEs with the algorithm of
Benettin et al. [20]. These vectors are usually called Gram-
SchmidtorbackwardLVs[21]andtheyspanthesubspacesin
tangentspacethat,atpresenttimet, havegrownatexponen- IV. LONG-RANGECORRELATIONSOFLYAPUNOV
tialratesλ sincetheremotepast. BackwardLVs{b (x,θ)} VECTORS
n n
dependontheparticularscalarproductadoptedfortheGram-
SchmidtorthogonalizationinBenettin’salgorithm(thoughthe Inthissectionwestudytheformofbothbackward(Gram-
subspaces they span are genuine). This and other undesired Schmidt)and characteristicLVscorrespondingto the largest
propertiesrenderthebackwardLVsunsuitedtoanalyzeprob- LEs, their localization properties, and the existence of long-
lemslikeextensivityorhyperbolicityquestions. rangecorrelations.
RuelleandEckmann[13,22]noticedtimeagothatonecan
LVsarestronglocalizedobjects,therefore,likeinRefs.[11,
defineadifferentsetofLVs,theso-calledcharacteristic[21]
12], it is very convenient to work with the associated fields
LVs{g (x,θ)},thatareintrinsictothedynamicsandareuni-
n obtainedafteralogarithmictransformation:
vocally defined (i.e. independent of how the scalar product
is defined). Each of these vectors is covariant with the dy-
h (x,θ)=ln|g (x,θ)| (5)
namics, and the corresponding n-th LE is indeed recovered n n
by either a forward or backward integration of an infinitesi-
malperturbationinitially alignedwith the n-thcharacteristic (andlikewise for backwardLVs). We will referto hn as the
LV. Actually, the perturbationwill be alignedwith the LV at nth surface due to the similarities of its dynamics with that
alltimes(hencethecovariance). HoweverforbackwardLVs of the kinetic rougheningof fractal surfaces in growth mod-
(other than the main one, because b = g ) the correspond- els(seedetailsbelow). Alsorecallthatthenormisirrelevant
1 1
ingLEsareonlyrecoveredintegratingbackward. Untilvery becauseaLVonlyindicatesadirectionintangentspace,and
recently, no efficient numerical algorithms were available to thisimpliesthatthe meanheightofhn is irrelevant(onlyits
computecharacteristicLVsinlargeextendedsystems. profile fluctuationsmatter); this shouldbe beard in mind for
In this paper we calculate the set of characteristic LVs by thetheoreticalanalysisbelow.
adapting to DDSs the method recently introduced by Wolfe Figure 1(a) shows a snapshot of the characteristic LV
and Samelson [14]. In order to find the nth characteristic surfaces h , h , and h corresponding to the three largest
1 2 3
3
0 backward LVs
hn−−−321000 t h1 h 2 h3 t+ τ/2 (a t)+ τ 106 ~ k−2
Time
104
− h1 100 h2 − h1 h3 − h1 ~ k−1
hn−10 t t+ τ/2 (b t)+ τ S(k)n 102 local slope −1
Time −1.2
−1.4
FIG. 1: (Color online) (a) Snapshot of the surfaces corresponding 100 −1.6
−1.8
to characteristic LVsn = 1,2,3 (τ = 3277 t.u.). (b) Difference −2
of2ndand3rdLV-surfaceswithrespecttothe1stLV-surface. The −2 −1
existenceof plateaus evidences thepiecewisecopy structureof LV 10−2 10 10
surfaceswithrespecttothe1stone. 10−3 10−2 10−1 100 101
k
characteristic LVs
LEs2.Wefindthatnon-leadingLVsurfacesareapproximately 106
piecewisecopiesoftheleadingLVsurfaceh1, seeFig.1(b). ~ k−2
This‘replicationproperty’isahighlynontrivialphenomenon
that was originally discovered to occur in chaotic extended 104
dissipativesystems[23,24]. Wealsoemphasizethatcharac- ~ k−1.2
teristicLVsexhibitatendencytoclusterize,contrarytoback-
wardLVswhoselocalizationsitesarescatteredduetotheim- (k) 102 local slope
Sn −1
posed orthogonality. Indeed, as exemplifiedby the snapshot −1.2
shown in Fig. 1(a), the first three characteristic LVs localize −1.4
(reachtheirlargestmagnitude)atthesamepoint(completely 100 −1.6
−1.8
overlappingontheleftmostregioninthisparticularsnapshot). −2
Thisdoesnotoccurallthetimebutinanintermittentmanner. −2 −1
10−2 10 10
Interestingly,thesestructuralfeatures–namely,replication
10−3 10−2 10−1 100 101
and clustering– have recently been reported to occur generi-
k
cally for LVs in chaotic spatially extendedsystems [23–25].
ThisdeepensintheanalogybetweenDDSsandsystemswith
FIG. 2: Structure factors of LV-surfaces for the MG model with
extensivechaosinonedimension,whichhappenstoholdeven
τ = 3277 t.u. The curves from top to bottom correspond to
atthelevelofnon-leadingLVs.
n = 1,4,8,16,32. The insets show the dependence of the local
Regardingtheexistenceoflong-rangecorrelations,wefind slopes of different curves. For n = 1 the asymptotic slope is −2
that LV surfaces have a self-affine spatial structure at long whereasforn>1theasymptoticslopesare≈−1and≈−1.2for
scales, which translates into power-law correlations. A de- backwardandcharacteristicLVs,respectively.
tailedanalysisofthespatialstructurecanbebestachievedby
computing the Fourier transform of the surfaces at discrete
times θ given by hˆ (k,θ) = 1 τ exp(ikx)h (x,θ)dx,
n √τ 0 n
with wavenumbers k ∈ [2π, π],Rwhere τ is the system
τ dt
size in this representation. For a given LV surface h ,
n called roughness exponentα = 1/2 (in agreement with the
the Fourier transform of the two-point correlator hh (x +
n 0 resultin[12] forthe Ikedamodel). Figure2 also revealsthe
x,θ)h (x ,θ)i−hh (x ,θ)i2,istheso-calledstructurefac-
n 0 n 0 existenceofaparticularcrossoverwavenumberk ≈ 2πn for
tor: n× τ
thenthLV.Theexponentsobservedfork <k are−1.2and
n×
S (k)=hhˆ (k,θ)hˆ (−k,θ)i (6) −1forcharacteristicandbackwardLVs, respectively. These
n n n
numberscoincidewiththosepreviouslyreportedindissipative
where the brackets denote average over time θ and realiza- systemsinonedimension[23,24]. Ourinterpretationofthese
tions. Therefore, the structure factor Sn(k) directly informs resultsisthatforbackwardLVsthek−1dependenceseemsto
aboutthenthLVsurfacecorrelationsatscale1/k. betheconsequenceofresidualcorrelationswith ageometric
InFig.2itcanbeseenthatthestructurefactorforthefirst origin in the orthogonalityof the basis. The −1.2 exponent
LV-surface decays asymptotically as k (2α+1), with the so- ofcharacteristicLVsindicatestheyconveymoreinformation
−
amongdistantpartsofthesystemthanbackwardLVs;which
canbeunderstoodasaconsequenceofthecovarianceofchar-
acteristicLVsastheyareconsistentwithbothpastandfuture
2Forτ =3277t.u.thereare168positiveLEs.Then0-thLEvanishes,with evolutions [characteristic LVs with n > 1 are “saddle solu-
n0≃0.0514τ+0.79. tions”of(2)forcedby(1),see[24]].
4
V. MAINLYAPUNOVVECTORANDTHE variablesθ =(t−x)/τ,sothat(2)isnowwrittenas:
UNIVERSALITYCLASSQUESTION
∂ δy(x,θ)=uδy(x,θ)+vδy(x,θ−1), (8)
x
Aswehaveseenintheprevioussection,thesurfaceh as-
1 where x ∈ [0,τ) is the ‘spatial’ position and θ ∈ Z is the
sociatedwiththemainLVexhibitsscale-invariantcorrelations
discretetime. Wehaveseenintheprevioussectionthatade-
in the k → 0 limit. This strongly suggests that this LV sur-
scriptionintermsofsurfacesinsteadofvectorsthemselvesis
face should belong to one of the universality classes of sur-
more appropriate due to the strong localization of the latter.
facegrowth. ThefactthattheLVsurfacesobeyscalinglaws
Moreover,to relate Eq. (8) with one of the stochastic partial
andthatsystemswithspatiotemporalchaoscouldbedivided
differential equations modeling surface growth, we wish to
intoafewuniversalityclasses,accordingtothescalingofthe
approximatethe discrete differencesin θ by a partialderiva-
associatedsurfaces,hasimplicationsinourunderstandingof
tive.Insum,wehavetotaketwosteps:
chaosin extendedsystemsfroma statisticalphysicspointof
view.Moreover,thegenericreplicationandclusteringofchar-
(i) Transformtoasurface:h(x,θ)=ln|δy(x,θ)|.
acteristicLVsalongthemainvectordirectioninextendedsys-
tems[23–25](alsoobservedherefortime-delaysystems,see
(ii) Approximateθ byacontinuousvariable: s(θ)−s(θ−
previoussection)makesitevenmoreinterestingtodetermine
1)→∂ s.
the universality class of the main LV, since this is expected θ
to provide a great deal of information about the structure of
Notethatthesetwostepsdonotcommute.
spaceandtimecorrelationsoftheLVcorrespondingtolead-
ingaswellassub-leadingunstabledirections.
For a wide class of one-dimensional spatio-temporal
1. LinearZhangmodel
chaoticsystems(see [12]) thestatistical featuresofthemain
LV surface are well captured by the 1d-Kardar-Parisi-Zhang
(KPZ)[27]stochasticsurfacegrowthequation: ThetreatmentofEq.(8)bySa´nchezet. al.[11]proceeded
with step (ii) before step (i). In more detail, after step (ii),
∂ h(x,t)=χ(x,t)+κ(∂ h)2+ν∂ h, (7) Eq.(8)becomes:
t x xx
where χ is a white noise. The KPZ equation defines itself ∂θδy =−(1/v)∂xδy+(u/v+1)δy, (9)
animportantuniversalityclassofsurfacegrowth. Regarding
chaoticsystems,inmanyofthemthemainLVsurfacebelongs andnowtransformingtothesurfacepicture[step(i)]onehas:
to the KPZ class. This includes coupled-map lattices, cou-
pledsymplecticmaps,thecomplexGinzburg-Landaumodel, ∂θh(x,θ)=−(1/v)∂xh+(u/v)+1, (10)
Lorenz 96 model, among others [12, 23, 24, 28]. Only LVs
sothatwegetanequationforthefieldh(x,θ). Thisequation
ofanharmonicHamiltonianlatticesareknown[29]toexhibit
wasalreadyanalyzedinRef.[11]andwesummarizesomeof
correlationsthatareclearlyinconsistentwithKPZexponents
itspropertiesin thefollowing. On theonehand,therandom
(α>α =1/2).
KPZ
drift term (1/v)∂ h gives rise to an effective diffusion term
Statistical mechanics teaches us that, given the fact that x
D∂ h. Ontheotherhand,thepresenceofv inthedenomi-
KPZrepresentsadynamicaluniversalityclass,onewouldex- xx
natorinduceslargefluctuationsofthosetermsproportionalto
pectthatalargecollectionofdifferentsystemscouldbelongto
ζ = 1/v every time that v takes values close to zero. Usu-
theKPZuniversalitydespitetheirapparentdifferences. Only
ally,theprobabilitydensityfunctionisalgebraicatthelowest
symmetries and conservationlaws would determine the uni-
order:
versality class. That may explain why Hamiltonian (energy
conserving)systems,incontrasttodissipativesystems,donot
P(v)∼|v|σ (|v|≪1), (11)
generallybelongtotheKPZuniversality.Onthisbasis,DDSs
were also proposed[11] to belongto a differentuniversality
with σ > −1 to be normalizable (and σ = 0 in gen-
class, namely the Zhang model, arguing that these systems
eral). In turn, P(ζ) is heavy tailed with large fluctuations
break the x → −x symmetry in the spatial representation.
P(ζ) ∼ |ζ| (2+σ) (thisis indeedin agreementwith datawe
However,thisconclusionwasincontradictionwithanearlier −
havecollectedinoursimulationsoftheMG model: σ ≈ 0).
work[12]whereKPZwaspostulated.Ouraiminthissection
ThesereasoningledtheauthorsofRef.[11]toconcludethat
is to clarify this question by means of a theoretical analysis
the main LV surface in DDSs behaves following (10) and
andextensivenumericalsimulations.
genericallyfallsintotheuniversalityclassofthelinearZhang
model[30]:
A. Theoreticalanalysis:KPZvs.Zhangequation
∂ h(x,θ)=D∂ h+c+ξ(x,θ), (12)
θ xx
Ourstartingpointistheevolutionequationintangentspace that describes surface growth driven by a heavy-tailed noise
(2), of which the main LV is the asymptoticsolution. In the ξ(x,θ) with a probabilitydistribution P(ξ) ∼ |ξ| (1+µ) for
−
spatialrepresentation(see Sec. II) we performthe changeof |ξ|≫1andtheindexµ=2(1+σ)(≥2ifσ ≥0).
5
2. KPZuniversality

τ=3277
Inthisworkwe proposeanalternativeanalysisofEq.(8), 2 τ=13107
10
where step (i) is carried out before step (ii). This would be τ=52429
moresuitedbecauseoneappliesthenonlineartransformation
(i)priortotheapproximation(ii).Thusafterdefiningthecor- θ)
responding LV surface, h(x,θ) = ln|δy(x,θ)|, Eq. (8) be- 2W(
comes 101
∂ h(x,θ)=u±vexp[h(x,θ−1)−h(x,θ)], (13)
x
wherethechoiceofsign±comesfromtheabsolutevalueand
is irrelevant for the arguments that follow. We approximate 100
the differenceh(x,θ −1)−h(x,θ) by thepartialderivative 100 101 102 103 104 105
−∂ h.Withnofurtherapproximations3andtakinglogarithms θ = t / τ
θ
inEq.(13)weget:
∂ h(x,θ)=ln|v|−ln|∂ h−u| (14)
θ x
2/3
where the term η = ln|v| is again a fluctuating noise-like 0.6
sourcedueto thechaoticcharacterof thetrajectories. Inthe e
v
neighborhoodofv =0theprobabilitydistributionofvis(11), ati 0.5
v
a−n∞d )th∼usel−ar(gσe+1v)a|ηlu|,eis.eo.fth|eη|noairseeeixspnoontehnetaiavlylytarialered:. P(η → c deri 0.4
We canobtaina moreintuitiveequationbymakinguseof mi 0.3
h
thesmallgradientapproximation,|∂xh|≪1.Expanding(14) arit 0.2
intheform: g
o
L
∂ h (∂ h)2 0.1
∂ h(x,θ)=η−ln|u|+ x + x +O[(∂ h)3], (15)
θ u 2u2 x
0
wheretheprototypicalquadraticterm(∂ h)2ofKPZappears. 100 101 102 103 104 105
x
Ifuisfluctuating,thedriftterm(∂ h)/uagainyields(atlarge θ = t / τ
x
scales, after a spatial averaging) an effective diffusion term.
Inthiscase,Eq.(15)wouldleadtotheKPZEq.(7). Incon- FIG. 3: (Color online) Evolution of the width from a random per-
trast, if u is a constant, like for instance in the case of the turbationforthreedifferentvaluesofτ. Thelowerpanelshowsthe
MGmodel,adiffusiontermwouldalsoeventuallyappearfar localderivative. Forthelargestvalueofthedelay, τ = 52429t.u.,
fromthesmallgradientlimitasaresultoftheabsolutevalue thelocalderivativeexhibitsaplateauabout0.63, closetothetheo-
reticalvalueforKPZ(2β = 2/3). Thecurvesareaveragesoverat
inln|∂ h−u|in(14). Finally,noticethatthenoise-liketerm
x
least1000realizations.
ln(|v|/|u|),whichisdifferentfromthatin(12),willnotlead
torareeventsanditsdistributionwillbeexponentiallydecay-
ingingeneral.
widthatagiventimewhenstartedfromaflatinitialperturba-
At variance with Sa´nchez et al. derivation [11]
tion [h(x,0) ≈ const.]. To do so we let randominitial per-
(cf. Sec. VA1) we have arrived here at the LV surface
turbationsδy(x,0)toevolveandmeasuretheaveragesurface
equation (15) by expanding on the local slope ∂ h, which
x widthattimeθas
is generally expected to be a reasonable approach when the
roughness exponent α < 1 since the spatial average over a
windowofextentℓshouldscaleash|∂ h|i ∼ℓα 1andgoes W2(θ)= h(x,θ)−h(θ) 2 , (16)
x ℓ − (cid:28) (cid:29)
tozeroaswecoarse-grainℓ≫1. (cid:2) (cid:3)
whereh(θ) = τ 1 τ h(x,θ)dxisthemeansurfaceposition
− 0
B. Growthexponent and brackets denotRe averaging over different initial realiza-
tions. The growth exponentβ is defined as the exponentof
the transient power-law growth: W2 ∼ θ2β, before the sur-
Inthefieldofgrowinginterfacesitiscustomarytoquantify
face fluctuations saturate (for θ ≪ θ ); while W2 ∼ τ2α,
thetemporalfeaturesbylookingatthegrowthofthesurface whensaturationisreached(forθ ≫θ×)andtheperturbation
isvirtuallyalignedwiththeLV:δy ∝×g ⇒h=h +const.
1 1
Thekeypointnowisthatβequals1/3forKPZ,whileitis1/4
fortheZhangmodelwithµ≥2[31].Thismakesthetimeex-
3Consideringthetimederivative∂ hsmallenoughsothattheexponential
θ ponentβ anexcellentindextodistinguishbetweenKPZand
inEq.(13)can beexpanded andonlythe lowest order mayberetained
Zhangbehavior(notethatthespatialexponentαisthesame
yields Eq. (10). Wepropose here a different derivation that avoids this
assumptionand,remarkably,leadstoadifferentresult. inbothmodels:α=1/2).
6
forallq. Fortrulyextremeeventdominatedfluctuations,this

lengthscaleisexpectedtodiverge(howeverslowly)withthe
101 system size, so that the fluctuation distribution is truly non-
Gaussian in the thermodynamiclimit. This slow divergence
was indeed measured for models in the Zhang universality
q=5
class[26].Inthecaseoftime-delaysystems,thiswouldcorre-
G(l)q 100 spondtohavinglc(τ)increasingwithτ. Incontrast,wefind
that the characteristic length l ∼ 10 remains constant even
c
τ=819 after an increase of the delay of 16 times. Again, the sim-
ulations by Sa´nchez et. al. [11] were carried out in systems
τ=3277
10−1 q=1 τ=13107 withdelaysthatweretooshorttoobtainconclusiveevidence
onthedependenceofl withτ. Certainly,wecannotruleout
c
10−1 100 101 102 103 104 for sure an extremely weak dependence of lc on τ, but if it
l (time units) exists it must be sub-logarithmicand well beyondthe preci-
sion that we can reach in our simulations. We concludethat
multiscalingoftheLVsurfacefluctuationsinDDSsseemsto
FIG. 4: (Color online) Height-height correlation functions G (l),
q
q = 1,...,5. Multiscalingisobservedatsmalll. Aboveacertain beashortscalephenomenonthathasnoeffectinthethermo-
characteristicsizel ∼ 10nomultiscalingisobserved. Notethatl dynamiclimit where the universalityclass is defined; in this
c c
isinsensitivetoτ,andthecurvesoverlapalmostperfectly(exceptat casedescribedbytheKPZequation.
largelduetofinite-sizeeffects).
VI. CONCLUSIONS
Figure 3 shows that W2 grows with an exponent that,
for large enough delays, progressively approaches the KPZ
In this work we have implemented (to our knowledge for
growth exponentβ = 1/3, which stronglysupportsthe
KPZ thefirsttime)characteristicLVsinDDSs. Adaptationofthe
KPZ asymptotic (τ → ∞) scaling for chaotic DDSs. In
method proposed in [14] to this kind of systems —together
the paper by Sa´nchez et al. [11] the largest time delay used
withthecomputercapabilitiesnowadaysavailable—allowed
τ = 4000t.u.resultedinanestimationβ = 1/4. Inthelight
us to reach systems with fairly largedelays, which serves to
ofoursimulationswithmuchlargerdelaysweconcludethat
investigatethe “thermodinamiclimit” of these systems. Our
thevalueofτ usedinRef.[11]wasinsufficienttodetectthe
resultsfortheLVscoincidequantitativelywiththoseobtained
trueasymptoticscalingexponent.
inextendeddissipativesystems[23,24].
Inadditionwehaverevisitedthequestionofwhichuniver-
sality class the main LV belongs to. After simulations with
C. Multiscaling
very large delays we may conclude that the main LV sur-
facefallsintotheuniversalityclassoftheKardar-Parisi-Zhang
Inadditiontothegrowthexponent,Ref.[11] alsoinvoked
equation. Our theoretical argumentssupport this conclusion
the existence of multiscaling of the main LV surface as an
aswell.
argument supporting its identification with the universality
Insum,DDSareequivalenttoextendeddynamicalsystems
classofthelinearZhangmodel(12). Thismodel,contraryto
in the sense that infinitesimal perturbationsexhibit the same
theKPZequation,leadstosurfacesthatexhibitmultiscaling,
exponents characterizing spatiotemporal correlations. DDS
i.e. strongly non-Gaussian tails, induced by extreme events.
have beentraditionallyconsideredto be differentbecause of
If we compute the qth-height-heightcorrelation function for thelackofextensivityoftheLyapunovspectrum:thepositive
pointsseparatedadistancel, exponentsapproachzeroas∼ 1/τ [7]andthe(Kolmogorov-
Sinai)entropysaturateswithτ. Howevertheidentificationof
1/q
G (l)= |h (x+l,θ)−h (x,θ)|q , (17) τ withasizeimpliesthatcomparisonsshouldbedoneintem-
q 1 1
D E poralunitsofθ = t/τ,andextensivityisthenrecovered. As
(wheretheoverlinedenotestheaverageoverxandthebrack- λnt = λnτθ = Λnθ, the redefined LEs Λn = τλn do not
etsdenotearealizationsaverage)onefindsthatGq(l) ∼ lαq, decay to zero with τ, and the Lyapunovspectrum converges
andmultiscalingexistsiftheroughnessexponentsα depend toadensityinthethermodynamiclimit.
q
ontheindexq.
We have computedthe qth orderheight-heightcorrelation
function G (l), given by Eq. (17), for several values of τ Acknowledgments
q
(i.e. system sizes); see Fig. 4. The distribution of the sur-
facefluctuationsshowsclearsignsofmultiscalingwithlocal D.P. acknowledges support by CSIC under the Junta de
roughnessexponentsα thatdependonqforlengthscalesbe- Ampliacio´n de Estudios Programme (JAE-Doc). Financial
q
lowatypicalscaleofaboutl ∼10.Thesenon-Gaussianfea- supportfromthe Ministeriode Cienciae Innovacio´n(Spain)
c
turesdisappearatlargerscales(l ≫ l ),wheretheα ≈ 1/2 underprojectNo.FIS2009-12964-C05-05isacknowledged.
c q
7
[1] T. Erneux, Applied Delay Differential Equations (Springer, [18] V.S.Udaltsov,J.-P.Goedgebuer,L.Larger,andW.T.Rhodes,
NewYork,2009). Phys.Rev.Lett.86,1892(2001).
[2] Delayed complex systems in Phil. Trans. R. Soc. A 368, [19] W.H.Press,S.A.Teukolsky,W.T.Vetterling,andB.P.Flan-
No. 1911, edited by W. Just, A. Pelster, M. Schanz, and E. nery, Numerical Recipes in Fortran 77: The Art of Scientific
Scho¨ll(2010). Computing, 2nded.(CambridgeUniversityPress,Cambridge,
[3] F.T.Arecchi,G.Giacomelli,A.Lapucci,andR.Meucci,Phys. 1992).
Rev.A45,R4225(1992). [20] G.Benettin,L.Galgani,A.Giorgilli,andJ.-M.Strelcyn,Mec-
[4] G.Giacomelli, R.Meucci, A.Politi, and F.T. Arecchi, Phys. canica15,9(1980).
Rev.Lett.73,1099(1994). [21] B. Legras and R. Vautard, in Proc. Seminar on Predictability
[5] G.GiacomelliandA.Politi,Phys.Rev.Lett.76,2686(1996). Vol.I,ECWFSeminar,editedbyT.Palmer(ECMWF,Reading,
[6] M. Wolfrum and S. Yanchuk, Phys. Rev. Lett. 96, 220201 UK,1996),pp.135–146.
(2006). [22] D.Ruelle,Publ.Math.IHES50,27(1979).
[7] J.D.Farmer,PhysicaD4,366(1982). [23] I. G. Szendro, D. Pazo´, M. A. Rodr´ıguez, and J. M. Lo´pez,
[8] M.LeBerreetal.,Phys.Rev.A35,4020(1987). Phys.Rev.E76,025202(R)(2007).
[9] B.Dorizzietal.,Phys.Rev.A35,328(1987). [24] D. Pazo´, I. G. Szendro, J. M. Lo´pez, and M. A. Rodr´ıguez,
[10] S.Lepri,G.Giacomelli,A.Politi,andF.T.Arecchi,PhysicaD Phys.Rev.E78,016209(2008).
70,235(1993). [25] M. Romero-Bastida, D. Pazo´, J. M. Lo´pez, and M. A.
[11] A. D. Sa´nchez, J. M. Lo´pez, M. A. Rodr´ıguez, and M. A. Rodr´ıguez,Phys.Rev.E82,036205(2010).
Mat´ıas,Phys.Rev.Lett.92,204101(2004). [26] A.-L.Baraba´sietal.,Phys.Rev.A45,R6951(1992).
[12] A.PikovskyandA.Politi,Nonlinearity11,1049(1998). [27] M.Kardar,G.Parisi,andY.-C.Zhang,Phys.Rev.Lett.56,889
[13] J.-P.EckmannandD.Ruelle,Rev.Mod.Phys.57,617(1985). (1986).
[14] C.L.WolfeandR.M.Samelson,Tellus59A,355(2007). [28] A.S.PikovskyandJ.Kurths,Phys.Rev.E49,898(1994).
[15] M.C.MackeyandL.Glass,Science197,287(1977). [29] A.PikovskyandA.Politi,Phys.Rev.E63,036207(2001).
[16] K.Ikeda,Opt.Commun.30,257(1979). [30] Y.-C.Zhang,J.Phys.(France)51,2129(1990).
[17] J.-P.Goedgebuer,L.Larger,andH.Porte,Phys.Rev.Lett.80, [31] C.-H.LamandL.M.Sander,Phys.Rev.E48,979(1993).
2249(1998).