26
On the ‘‘genetic memory’’ of chaotic attractor of the barotropic ocean model V. Dymnikov a , Ch. Kazantsev b , E. Kazantsev c, * a Institute of Numerical Mathematics, 8, Gubkin str., 117333 Moscow, Russia b D epartement de Math ematiques, Universit e de Nancy 1, BP 239, 54506 Vandoeuvre-l es-Nancy, France c Projet Numath, INRIA Lorraine, 615 rue du Jardin Botanique, BP101, 54600 Villers-l es-Nancy, France Accepted 17 September 1997 Abstract The structure of attractor of barotropic ocean model is studied in this paper. Theorems of the existence of the attractor for the finite dimensional approximation of this model are proved as well as its convergence to the attractor of the model itself. Some properties of stationary solutions of this model and their stability are discussed. The structure of the attractor is partially explained by the sequence of bifurcations the system is subjected to by variations of leading parameters. The principal feature of the studied system is the existence of two ‘‘almost invariant’’ basins of chaotic attractor with very rare transitions between them. This is related to the rise of a couple of non-symmetric stable stationary solutions in the model with symmetric forcing. The ‘‘memory’’ of chaos appears also in the presence of maxima in the spectrum of energy. These maxima correspond either to the principal frequency of the limit cycle which arose in the Hopf bifurcation, or to the frequencies of the Feigenbaum pheno- menon. Ó 2000 Published by Elsevier Science Ltd. All rights reserved. 1. Introduction The study of the low frequency variability of the atmospheric and oceanic circulations is the question of principal importance because it provides us with the scientific basis for understanding and estimating of their predictability. If real climatic system is supposed to be described by a nonlinear dissipative model with a global attractor [1,2], then the study of the attractor structure, the dynamics of the system on the attractor and the stability of the attractor structure are problems of key importance in the theory of climate. As it has been shown in recent papers [3–5], even simple models of general circulation of the atmosphere and ocean exhibit the phenomenon of dynamical chaos on the attractor. This phenomenon occurs due to local instability of trajectories. One can use global Lyapunov exponents, Kolmogorov entropy or attractor dimension as global at- tractor characteristics. In spite of the fact that these characteristics may help to understand the dynamics and predictability of the system, they must be considered as a part of the problem only. The second part concerns the study of the probabilistic measure concentrated on the attractor of the climatic system. Non- uniformity of the measure on the attractor, or in other words, the presence of ‘‘preferable states’’, or clusters on the attractors allows us to classify the attractor’s structure as from the point of view of the predictability of trajectories in clusters, and from the practical point of view for real climatic systems. Such a system has been a subject of studies in meteorology and oceanography for a long time. www.elsevier.nl/locate/chaos Chaos, Solitons and Fractals 11 (2000) 507–532 * Corresponding author. 0960-0779/00/$ - see front matter Ó 2000 Published by Elsevier Science Ltd. All rights reserved. PII: S 0 9 6 0 - 0 7 7 9 ( 9 9 ) 0 0 1 0 2 - 2

On the “genetic memory” of chaotic attractor of the barotropic ocean model

Embed Size (px)

Citation preview

Page 1: On the “genetic memory” of chaotic attractor of the barotropic ocean model

On the ``genetic memory'' of chaotic attractor of the barotropicocean model

V. Dymnikov a, Ch. Kazantsev b, E. Kazantsev c,*

a Institute of Numerical Mathematics, 8, Gubkin str., 117333 Moscow, Russiab D�epartement de Math�ematiques, Universit�e de Nancy 1, BP 239, 54506 Vandoeuvre-l�es-Nancy, France

c Projet Numath, INRIA Lorraine, 615 rue du Jardin Botanique, BP101, 54600 Villers-l�es-Nancy, France

Accepted 17 September 1997

Abstract

The structure of attractor of barotropic ocean model is studied in this paper. Theorems of the existence of the attractor for the ®nite

dimensional approximation of this model are proved as well as its convergence to the attractor of the model itself. Some properties of

stationary solutions of this model and their stability are discussed.

The structure of the attractor is partially explained by the sequence of bifurcations the system is subjected to by variations of leading

parameters. The principal feature of the studied system is the existence of two ``almost invariant'' basins of chaotic attractor with very

rare transitions between them. This is related to the rise of a couple of non-symmetric stable stationary solutions in the model with

symmetric forcing.

The ``memory'' of chaos appears also in the presence of maxima in the spectrum of energy. These maxima correspond either to

the principal frequency of the limit cycle which arose in the Hopf bifurcation, or to the frequencies of the Feigenbaum pheno-

menon. Ó 2000 Published by Elsevier Science Ltd. All rights reserved.

1. Introduction

The study of the low frequency variability of the atmospheric and oceanic circulations is the question ofprincipal importance because it provides us with the scienti®c basis for understanding and estimating oftheir predictability. If real climatic system is supposed to be described by a nonlinear dissipative model witha global attractor [1,2], then the study of the attractor structure, the dynamics of the system on the attractorand the stability of the attractor structure are problems of key importance in the theory of climate.

As it has been shown in recent papers [3±5], even simple models of general circulation of the atmosphereand ocean exhibit the phenomenon of dynamical chaos on the attractor. This phenomenon occurs due tolocal instability of trajectories.

One can use global Lyapunov exponents, Kolmogorov entropy or attractor dimension as global at-tractor characteristics. In spite of the fact that these characteristics may help to understand the dynamicsand predictability of the system, they must be considered as a part of the problem only. The second partconcerns the study of the probabilistic measure concentrated on the attractor of the climatic system. Non-uniformity of the measure on the attractor, or in other words, the presence of ``preferable states'', orclusters on the attractors allows us to classify the attractor's structure as from the point of view of thepredictability of trajectories in clusters, and from the practical point of view for real climatic systems. Sucha system has been a subject of studies in meteorology and oceanography for a long time.

www.elsevier.nl/locate/chaos

Chaos, Solitons and Fractals 11 (2000) 507±532

* Corresponding author.

0960-0779/00/$ - see front matter Ó 2000 Published by Elsevier Science Ltd. All rights reserved.

PII: S 0 9 6 0 - 0 7 7 9 ( 9 9 ) 0 0 1 0 2 - 2

Page 2: On the “genetic memory” of chaotic attractor of the barotropic ocean model

The fact of existence of ``preferable'' clusters on the chaotic attractor or maxima in power spectra is wellknown (see for example [26,27]). However, the nature of these formations is not clear. Based on the simplebarotropic ocean model in a square, we try to answer the question ± does the attractor with chaotic dy-namics ``remember'' the chain of bifurcations it is generated by. This ``memory'' can help to associateclusters and maxima in power spectrum with stable states which disappeared during bifurcations.

This question is not completely new, one can ®nd remarks on this subject in [6±8], however, no paperproves one deals with a chaotic attractor of the system.

In this paper we consider a system for which one can prove that there exists a stable equilibrium onlywhen the external forcing is su�ciently small. This equilibrium forms the global attractor of the system[9,10]. When the forcing becomes larger, there can appear multiple equilibria. The number of equilibriamust be odd [11]. Due to possible Hopf bifurcation, they can loose their stability and generate limit cycles.However, one can suppose these equilibria leave some traces on the global attractor of the system. One cansuppose also that if the transition to chaos is realized by the Feigenbaum phenomenon of doubling period,the frequency of the original limit cycle or the multiple ones can be identi®ed in the ®nal power spectrum.

Although the model we use is a quite simple model of the barotropic ocean in a square with a symmetricdistribution of the wind stress, we can suppose the results can relate to the problem of multiple regimes ofboundary currents in the real ocean [6]. Our results are close to the results obtained in the paper [6],however our main purpose is the study of global attractor with chaotic dynamics. We should note also, thelow resolution used for the model allows us to interpret the results as related to some general case thanks tothe fact they are close to former results with high resolution models [6] and thanks to the phenomenon of``spurious chaos'' [4] which can be presented in this model.

The paper is organised in the following way. In the ®rst section, the di�erential model of the barotropicoceanic circulation is described. Its ®nite dimensional approximation, method of time integration andcalculation of Lyapunov exponents and attractor dimension are presented in Section 2. Section 3 is devotedto the proof of the theorem of existence of the global attractor of the model dynamics. Calculations of twobifurcation chains from stationary point to chaotic attractor for two di�erent model viscosities and dis-cussion are presented in Section 4.

2. Model equations and method of solution

2.1. Barotropic vorticity equation

In this paper we consider ocean dynamics in the barotropic formulation, i.e. all the thermodynamice�ects are neglected and the vertical structure of the ocean is considered as uniform. The equation of thedynamics of the wind-driven ocean is written in the usual way

oDwot� J�w;Dw� `� � lD2wÿ rDwÿ f ; �1�

where w is the barotropic streamfunction. Associated velocity ®elds can be calculated as

u � ÿ owox; v � ow

oy: �2�

Operators presented in the equation are Jacobian and Laplacian ones:

J�w;x� � owox

oxoy

�ÿ ow

oyoxox

�; Dw � o2w

ox2� o2w

oy2:

We assume the b-plane approximation for the Coriolis parameter f , which represents the e�ect of the Earthrotation in this equation, i.e. we suppose that this parameter is linear in the y coordinate:

` � `0 � b�y ÿ y0�; �3�where `0 and y0 are the value of the Coriolis parameter and the coordinate of the mid-latitude of the basin.We use very simple basin geometry represented by a square box of side length L � 4000 km. We supposethat this basin is located in the middle of the North Atlantic, so we take the value of the Coriolis parameter

508 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 3: On the “genetic memory” of chaotic attractor of the barotropic ocean model

in the middle of the basin to be equal to `0 � 9:3� 10ÿ5 sÿ1, and its meridional gradient b � 2�10ÿ11 �m� s�ÿ1

.The source of energy in this equation is presented by the atmospheric wind stress applied to the surface.

In this paper, we take a steady zonal wind with a now classical two gyre antisymmetric pattern [12]. This isseen as a schematic pattern for the mean curl of the wind stress over the North Atlantic ocean in middlelatitudes. Its magnitude is equal to

f � 2ps0

qHLsin

2pyL; �4�

where s0 is the wind tension on the surface, q � 1000 kg mÿ3 is the density of water. The depth of the activewind driven current H has been chosen as 500 m. The dissipation in the Eq. (1) is composed by the har-monic lateral friction AD2w and the bottom drag rDw.

The Eq. (1) is subjected to impermeability and free-slip boundary conditions

wjoX � 0; DwjoX � 0: �5�

2.2. Numerical approximation

In order to look for a weak solution of the problem (1) and (5), we perform its variational formulation.We denote the L2 scalar product by h:; :i:

hw;ui �Z Z

Xwudxdy: �6�

Let us multiply Eq. (1) by a test function u�x; y� 2 H 10 �X� and perform the integration over the domain X.

Using notation (6) we get

oxot;u

� �� hJ�w;Dw� by�;ui � ÿlhrx;rui ÿ rhx;ui ÿ hf ;ui;

hx;ui � ÿhrw;rui:�7�

Here we have used the Green formula:

hDw;ui � ÿhrw;rui �Z

oX

owon

uds � ÿhrw;rui

thanks to the fact that / 2 H 10 �X�, i.e. ujoX � 0.

The variational formulation (7) and (8) of the problem (1) allows us to look for a solution by the ®niteelement method (FEM). The feasibility and utility of FEMs for modelling ocean dynamics was ®rst ad-dressed in [13]. He stated that using FEMs brought numerous advantages of modelling like precision,natural conservation of model invariants, ¯exibility of discretisation of complex domains, etc. These ad-vantages remain even when irregular discretisation of the domain is performed.

So far, the solution produced by the barotropic model of the North Atlantic typically includes a westernboundary layer with intense velocity gradients; the advantage of re®ning the triangulation along the westernboundary of the domain is rather clear. This helps to keep the quality of explicit eddy resolution by themodel while working with a lower number of grid nodes. The comparison of ®nite element and ®nitedi�erence models performed in [14] revealed that the di�erence arising between simulations by FE and FDtechniques can be regarded as insigni®cant when the number of FE nodes is about six times lower than thenumber of FD ones.

In spite of the fact that the number of operations per time step and grid node is much higher for the FEmodel, the possibility of reducing the number of grid points considerably diminishes the computational costof a model run. The possibility of a better precision working with fewer grid points is very valuable in thiswork (especially for Lyapunov exponents calculations) due to a high number of operations per point.

The package MODULEF [15] has been used to perform a triangulation of a domain. This packageproduce quasi-regular triangulation of the domain based upon the prescribed grid nodes on its boundary.

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 509

Page 4: On the “genetic memory” of chaotic attractor of the barotropic ocean model

The grid used in this paper is presented in Fig. 1. It is the result of the triangulation of a square of unit sidelength. We require the re®ning of the triangulation near the western boundary, especially in the middle ofthe domain where velocity gradients are extremely sharps. This triangulation is composed of 92 triangles.The integration points set, being a union of vertices and centers of edges of triangles, counts 211 nodes. Sothe resolution of the grid varies between 1=40 of the side length (about 100 km) near the western boundaryand 1=7 of the side length (about 550 km) near the eastern one.

Following [16] we use P2 ®nite elements, i.e. polynomials of second degree pi�x; y� � aix2 � bixy�ciy2 � dix� eiy � fi. Finite elements pi�x; y� are constructed so that being of the second order they are equalto 1 at the ith integration point and zero at all other points.

According to the Dirichlet boundary conditions, we consider internal points of the domain only:�xi; yi� 2 X n oX for i � 1; . . . ;N :, so the functions w;x are presented as linear combinations

w�x; y; t� �XN

i�1

wi�t�pi�x; y�; x�x; y; t� �XN

i�1

xi�t�pi�x; y�: �8�

Using these expressions, the discretised system (7) and (8) is written:Xi

oxi

othpi; pji �

Xi

Xm

wm�xi � byi�hJ�pm; pi�; pji � ÿlX

i

xihrpi;rpji ÿ rX

i

xihpi; pji

ÿX

i

fihpi; pji; �9�Xi

xihpi; pji � ÿX

i

wihrpi;rpji: �10�

To simplify notations we de®ne matrices of mass and rigidity as

Mj;i � hpi; pji; Cj;i � hrpi;rpji i � 1; . . . ;N ;j � 1; . . . ;N :

��11�

Then the system (9) and (10) writes in matricial form

Moxot�X

i

Xm

wm�xi � byi�hJ�pm; pi�; pji � ÿlCxÿ rMxÿMf ;

Mx � Cw:

�12�

Fig. 1. Triangulation of a unit square.

510 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 5: On the “genetic memory” of chaotic attractor of the barotropic ocean model

This system has been forwarded in time as follows

Mxn�1 ÿ xnÿ1

2s�X

i

Xm

wnm�xn

i � byi�hJ�pm; pi�; pji

� ÿlCxn�1 � xnÿ1

2ÿ rM

xn�1 � xnÿ1

2ÿMf : �13�

The scheme (13) is of second order accuracy in time.

2.3. Lyapunov exponents

To verify the solution is chaotic with a strange attractor we calculate its Lyapunov exponents. They areintroduced to describe the growth or decay of in®nitesimal elementary m-dimensional volume around thepoint �wt�0�x; y� in the limit of large time.

ki � limT!1

1

Tln lim

e!0supF2H ;

dim F�i

infx2F ;kxt�0k<e

kxt�Tkkxt�0k : �14�

In®nitesimal limit norm of initial perturbation allows us to use linear approach. The operator of (12)linearised with respect to the trajectory �w�x; y; t� can be written

Moxot�X

i

Xm

�wm�D �wi � byi� � �wmxi�hJ�pm; pi�; pji � ÿlCxÿ rMx

Mx � Cw

�15�

with Dirichlet boundary conditions. We use Crank±Nicholson scheme for time integration

Mxn�1 ÿ xn

s�X

i

Xm

�wn�1=2m

xn�1i � xn

i

2

�� � �xn�1=2

i � byi�wn�1m � wn

m

2

�hJ�pm; pi�; pji

� �ÿlCÿ rM�xn�1 � xn

2�16�

or in the resolved form

xn�1 � B� �wn�1=2�xn:

After L time steps we obtain

xn�L � B� �wn�Lÿ1�xn�Lÿ1 � B� �wn�Lÿ1�B� �wn�Lÿ2�xn�Lÿ2 �YLÿ1

i�0

B� �wn�i�xn � B�L�� �wn�xn: �17�

It has to be noted here, that the index �L� is not the power of the operator. This index denotes the numberof multiplications of di�erent operators B� �wn�i�. This index corresponds to the length of the basic trajectory�w or to the time Ls during which we let the initial error to grow. Further, we shall address this time as thecomposition length.

To estimate the evolution of the norm of xt�T in (14) we de®ne the scalar product in the ®nite dimen-sional space of discretised functions. We shall use the ®nite dimensional approximation of the scalarproduct in L2 functional space:

hh;/i �Z

Xh/dxdy �

XN

i

XN

j

hi/j

ZX

pipj dxdy �XN

j

/j�Mh�j � �Mh;/�: �18�

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 511

Page 6: On the “genetic memory” of chaotic attractor of the barotropic ocean model

The norm of a function de®ned by this scalar product is

khk � hh; hi1=2 � �Mh; h�1=2:

Since the mass matrix M is positive de®nite, the scalar product �Mh; h� is positive for all h.We can obtain the norm of the perturbation after L time steps from (17):

kHt�Lsk2 � �MHt�Ls;Ht�Ls� � �MB�L�Ht�0;B�L�Ht�0� � �B�L�tMB�L�Ht�0;Ht�0�; �19�where B�L�

tis the transposed operator.

For our system the evolution of an error is expressed by (17) in discrete time. The time interval is equal tothe product of number of steps and the time step value T � Ls. The evolution of the norm of error is givenby (19), so we can calculate Lyapunov exponents as

ki � limL!1

1

2Lsln lim

e!0supF2H ;

dim F�i

infH2F ;kHt�0k<e

�B�L�tMB�L�Ht�0;Ht�0��MHt�0;Ht�0�

� supF2H ;

dim F�i

infH2F ;kHt�0k<e

�MMÿ1B�L�tMB�L�Ht�0;Ht�0�

�MHt�0;Ht�0� : �20�

This ratio is equal to the ith eigenvalue of the matrix Mÿ1B�L�tMB�L�. Therefore, we can calculate the

Lyapunov exponents ki as

ki � limL!1

1

2Lsln li�Mÿ1B�L�

tMB�L��; �21�

where li�Mÿ1B�L�tMB�L�� are eigenvalues of the problem

Mÿ1B�L�tMB�L�/i � li/i �22�

ordered as l1 P l2 P � � � P lN :Multiplicative ergodic theorem proved in [17] reveals the existence of the limit matrix

limL!1�Mÿ1B�L�tMB�L��1=2Ls

. More to the point, if the system is ergodic, Lyapunov exponents (21) do notdepend on trajectory (i.e. initial data, except the set of zero measure), they characterise the whole attractorof the system. In other words, if the global attractor is the only invariant set in the phase space, theLyapunov exponents characterise the attractor in which the trajectory lies rather than the trajectory itself. Ifthe phase space of the system is split into several local attractors, then ki can be quite di�erent in these localattractors. However, even in this case their values relate to characteristics of the local attractor de®ned bytrajectory.

There are some pieces of information we can extract from these exponents. First, if none is positive, thebehaviour of the system is not chaotic. In this case, a small perturbation to a trajectory remains near to orapproaches that trajectory. The sum of all Lyapunov exponents of any dissipative system is negative. If thissum is positive, the system will ``blow up'' in time. This fact indicates us that, probably, we apply wrongphysics or we have some kind of numerical instability. If this sum is zero, we have a conservative Ham-iltonian system and certainly no strange attractor.

Calculated Lyapunov exponents allow us to calculate the attractor dimension de®ned by the trajectory.Let J be the number of Lyapunov exponents which produce a positive sum, but adding kJ�1 we obtain thesum negative:

PJi�1 ki P 0; but

PJ�1i�1 ki < 0. Then the value of the attractor dimension we can ®nd as [18]

D � J �PJ

i�1 ki

jkJ�1j �23�

However, these exponents, which we shall call global Lyapunov exponents, relate to the in®nite time scalesof a system behaviour. In practice we are frequently interested in how perturbations grow or decay during®nite time, i.e in a ®nite number of time step L. In meteorological or oceanographical application it is muchmore important to estimate the predictability for a certain period of evolution, say one month, rather thanfor thousands years. So, in this paper we consider also the generalisation of the notion of Lyapunov

512 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 7: On the “genetic memory” of chaotic attractor of the barotropic ocean model

exponents introduced in [19]. Thus, we are able to address local growth rate of instabilities across theattractor. To study the growth or decay of perturbations over a ®nite lag time from any given point on theattractor we calculate local Lyapunov exponents of the composition length Ls:

k�L�i � �w� �1

2Lsln li�Mÿ1B�L�

tMB�L��: �24�

These exponents characterise the growth of a perturbation made at point �w after a ®nite number L timesteps evolution by the dynamics. If local exponents are large then the possibility to predict is rather re-stricted in this part of the phase space. On the other hand, if they are small or even negative (even whenglobal exponents might be positive), then the predictability in that part of the phase space is enhanced withrespect to the average one. Local Lyapunov exponents k�L�i � �w� depend on the particular trajectory �w (or itsinitial data) and on the period Ls during which the integration is performed.

We shall consider also averages of these exponents over attractor. As the approximation of this averagewe use averaging over large number of trajectories �wm started from di�erent initial points on the attractor.To be sure that initial points belong to the attractor we choose them as di�erent points on the same tra-jectory obtained by long time integration of model equations.

�k�L�i � limM!1

1

M

XM

m

k�L�i � �wm�: �25�

Like global exponents de®ned by (21), the �k�L�i do not depend on the trajectory, so they are characteristics ofthe attractor and not of some speci®c orbit on it. They show us how much, on the average over the at-tractor, a perturbation behaves in time Ls.

In [3] we have shown that when the time Ls grows· the largest average local Lyapunov exponent �k�L�i does not grow,· the average local Kolmogorov entropy �K�L�� �w� does not grow,· the total sum of average local Lyapunov exponents does not depend on Ls, i.e. remains constant.

Here we should note that numerical solution of eigenvalues problem (22) may not be easy due to the factthat matrix B�L�, being a product of L matrices, is very ill conditioned for high L. The condition numbercalculated for the spectral norm of a selfadjoint positive de®nite matrix is equal to the relation of the largesteigenvalue to the smallest one. So in our case this number is equal to

cond�Mÿ1B�L�tMB�L�� � l1

lN� exp2Ls�k1ÿkN � :

Thus, the condition number of this matrix depends exponentially on the space resolution scale and on thelength of the time interval Ls and can reach rather high values.

This leads to large errors in eigenvalues obtained by standard algorithms. To be capable to calculatelocal Lyapunov exponents with long Ls and to study the convergence of local exponents to global oneswhen Ls!1 we use the improved Eckmann and Ruelle method described in [19].

To verify that the process converges to the global Lyapunov exponents, we calculate averaged localexponents (25) for di�erent composition lengths Ls.

To increase the accuracy of estimates of global exponents, local values have been extrapolated to in®nity.The extrapolation procedure is based on the fact that average local exponents decrease linearly in loga-rithmic coordinates for su�ciently large L [19,3,20]. That means there exists a slanted asymptote in the�log L! �1� limit.

ln�k�L�i ÿ k1i � � Ai log2 L� Bi when log L! �1: �26�However, it is clear there exists another limit of Lyapunov exponents, when t! 0. This limit has beencalled instantaneous Lyapunov exponents k�I� in [21]. This limit can be represented by log L! ÿ1 in thelogarithmic coordinates.

ln�k�L�i ÿ k1i � � k�I�i when log L! ÿ1: �27�Hence, the curve ln�k�L�i ÿ k1i � possesses two asymptotes: slanted (26) in the �1 limit and a horizontal one(27) in the ÿ1 one. We suppose that this curve can be approximated by a hyperbole

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 513

Page 8: On the “genetic memory” of chaotic attractor of the barotropic ocean model

ln�k�L�i

�ÿ k1i � ÿ ln�k�I�i �

�Ai log2 L�

� Bi ÿ ln�k�L�i ÿ k1i ��� Const

and we use only three lower terms of the Taylor series

ln�k�L�i ÿ k1i � � Ai log2 L� Bi � Ci

log2 L: �28�

The coe�cients k1i ;Ai;Bi;Ci have been calculated by RMS method from average local exponents. Obtainedvalue of k1i is used as the extrapolation of the global exponent.

This extrapolation has been tested on the stable periodic solution, when the largest global Lyapunovexponent is known to be equal to 0. The extrapolated value (28) does not exceed �10ÿ7, while the ex-trapolation made by (26) provides the accuracy of �10ÿ3 only.

This extrapolation will be used in this paper to calculate the dimension of the attractor by the KaplanYorke formula (23).

3. Existence of global attractor

3.1. Stationary solutions of the QG model

Let us recall a general result about stationary solutions of ordinary di�erential equations. We consider asystem of equations [22]

_x� F �x� � 0; x 2 Rm; �29�where F is a su�ciently regular function such that �F �x�; x� > 0 for jxjP r and F �x� 6� 0 for jxj � r. Then wehave the following theorem:

Theorem 3.1. Let F �x� be a continuous function defined in Rm and such that for r sufficiently large, it verifies�F �x�; x� > 0 for jxjP r and F �x� 6� 0 for jxj � r. Then· The Eq. (29) admits at least one stationary solution,· The number of stationary solutions, such that the determinant of its Jacobian matrix det�oFi=oxj� is non

zero, is odd: N � 2M � 1; M P 0,· Among all N � 2M � 1 stationary solutions, at least M solutions are instable.

We consider the space Vm spanned by the m ®rst eigenfunctions uk; k � 1; . . . ;m of the Laplace operatorÿD with the Dirichlet boundary conditions. We look for a w 2 Vm that satis®es the Eq. (1) for f 2 Vm. Thedevelopment of w and f on the basis of Vm has a form w�t; x� �Pm

k�1 wk�t�uk�x�; f �x� �Pmk�1 fkuk�x�;

where uk�x� are normalised byR

X jrukj2 dx � 1. Then Dw�t; x� � ÿPmk�1 kkwk�t�uk�x�. The Eq. (1) in this

space can be written in form

ÿ dwk

dtÿ lkkwk ÿ rwk � b

Xm

i�1

wi

ZX

oui

oxuk dx�

Xm

i;j�1

wiwjkj

ZX

J�ui;uj�uk dx � fk

kk: �30�

We de®ne the vector F as

Fk � �lkk � r�wk ÿ bXm

i�1

wi

ZX

oui

oxuk dxÿ

Xm

i;j�1

wiwjkj

ZX

J�ui;uj�uk dx� fk

kk: �31�

The scalar product of w with the vector F can be written

F �w�;w� � �Xm

k�1

Fk�w�wk

�Xm

k�1

�lkb�Xm

i;k�1

wkwi

ZX

oui

oxuk dxÿ

Xm

i;j;k�1

wiwjwkkj

ZX

J�ui;uj�uk dx�Xm

k�1

wkfk

kk: �32�

514 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 9: On the “genetic memory” of chaotic attractor of the barotropic ocean model

A simple computation shows that

bXm

i;k�1

wkwi

ZX

oui

oxuk dx � b

ZX

owox

wdx � 0;

and Xm

i;j;k�1

wiwjwkkj

ZX

J�ui;uj�uk dx �Z

XJ�w;Dw�wdx � 0:

Then

F �w�;w� � �Xm

k�1

�lkk � r�w2k �

Xm

k�1

wkfk

kkPXm

k�1

�lkk � r�w2k ÿ

Xm

k�1

jwkjjfkjkk

PXm

k�1

�lkk � r�w2k ÿ

Xm

k�1

lkkjwkj2 ÿXm

k�1

jfkjlk3

k

�Xm

k�1

rw2k ÿ

Xm

k�1

jfkjlk3

k

So F �w�;w� � > 0 for

jwjP����������������������������1

r

Xm

k�1

jfkjlk3

k

� 1

s� R:

Hence, we apply the Theorem 3.1 and get

Theorem 3.2.

· The Eq. (1) admits at least one stationary solution,· The number of stationary solutions, such that the determinant of its Jacobian matrix det�oFi=oxj� is non

zero, is odd: N � 2M � 1; M P 0,· Among all N � 2M � 1 stationary solutions, at least M solutions are instable.

3.2. Attractor

We consider the system (1) with boundary conditions (5) for f 2 Hÿ1�X� and x0 2 Hÿ1�X�.oxot� J�w;x� `� � lDxÿ rxÿ f in X� R�;

x � Dw in X� R�;

x�x; 0� � x0�x� in X;

x � w � 0 on oX:

�33�

We equip Hÿ1�X� with the scalar product

hhx1;x2ii � ÿÿ1hx1;Dÿ1x2i�1 �

ZXrw1rw2 dx: �34�

It de®nes the norm k � kÿ1 which is equivalent to the usual norm in Hÿ1�X�. We denote by j � j and ��; �� theusual norm and scalar product in L2�X�. We recall the following theorems [9].

Theorem 3.3.

· The problem (33) admits a unique solution x�t; x� in C��0; T �;Hÿ1�X�� \ L2�0; T ; L2�X�� \ L2loc�0; T ;H 1

0 �X��for all T > 0.

· This solution verifies for t P 0� l

R T0jxj2 dt6K1T � jrw�0�j2

� jrw�t�j26 jrw�0�j2 exp�ÿ�r� lk1�t� � K22 1ÿ exp�ÿ�r� lk1�t�� �

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 515

Page 10: On the “genetic memory” of chaotic attractor of the barotropic ocean model

� for 0 < t < T jx�t�j26K3jrw�0�j2 � K4 where K1; K2; K3; K4 depend on X; l; r; f ; T .· The mapping x0 ! x�t; x� is continuous in Hÿ1�X�.

Theorem 3.4. The semi-group G�t� from Hÿ1�X� to Hÿ1�X�, associated to the problem (33), is such that· there exists a bounded absorbing set in Hÿ1�X�,· the semigroup S�t� associated to the system (33) is �Hÿ1�X�; L2�X��-bounded in Hÿ1�X�,· there exists a maximal attractor A which is bounded in L2�X�, compact and connected in Hÿ1�X�. Its basin

of attraction is the whole space Hÿ1�X�. This attractor has finite Hausdorff and fractal dimensions,· for any t0 > 0 S�t0� is uniformly differentiable. Its derivative at x0 is a linear operator L�t0;x0� in Hÿ1�X�.

Moreover, we have supx02A jL�t0;x0�jHÿ1 <1.

Let us now consider the space Vm, spanned by the ®rst m eigenfunctions of the Laplace operator inHÿ1

0 �X�. We denote the orthogonal projector of Hÿ1�X� onto Vm by Pm. The projection of the system (33)onto the space Vm has a form

dxm

dt� PmJ�w;x� `� � lDxm ÿ rxm ÿ Pmf ;

xm � Dwm; �35�xm�x; 0� � xm0�x� � Pmx0:

We introduce the map Sm�t� : Vm ! Vm de®ned by the semigroup Sm�t� associated with the system (35). Allthe results of the Theorem 3.4 hold for the system (35) [9]. Finally we recall some de®nitions and a theoremof [23].

De®nition 3.1. The semigroup Sh�t� conditionally approximates S�t� on U � Hÿ10 �X� uniformly on a

compact set I � �t0; t1� � R� if there are a constant h�I ;U� and a function g�h; I ;U� de®ned for0 < h6 h�I ;U� such that limh!0 g�h; I ;U� � 0 and for any 0 < h6 h�I ;U� if u 2 U \ Vh has the propertythat S�t�u; Sh�t�u are de®ned and belong to U for t 2 �0; t2�; t0 < t26 t1 then kS�t�uÿ Sh�t�ukÿ16 g�h; I ;U�for t0 < t6 t2.

De®nition 3.2. A semigroup S�t�; t P 0 in a Banach space is said to be asymptotically smooth if for anybounded set B � V there is a compact set J � J�B� � V such that J attracts the setfx 2 B; S�t�x 2 B for t P 0g. We denote by N�B; e� the e-neighbourhood of the set B in a Banach space.

Theorem 3.5. Assume that S�t� has a local compact attractor A and N1 denotes an open neighbourhood of Asuch that A attracts N1. Suppose that there are constants h0 > 0; d0 > 0; t0 P 0 and two open neighbourhoodsN2; N3 of A with N1 � N2 � N�N2; d0� such that for 0 < h6 h0,· S�t�N1 � N2 for t P 0,· Sh�t��N1 \ Xh� � N2 for 0 < t6 t0,· for any xh 2 N�N2; d0� \ Vh, there exists T �xh� > 0 such that Sh�t�xh 2 N3 for 0 < t6 T �xh�,· Sh conditionally approximates S�t� on N3 uniformly on compact set of �t0;1�Then, if each Sh�t� is asymptotically smooth, there is h0 > 0 such that, for 0 < h�t�6 h0, Sh�t� admits a localcompact attractor Ah, which attracts N1 \ Xh. Moreover, dV �Ah;A� ! 0 as h! 0, where dV �A;B� �supx2A distV �X ;B� � supx2A inf y2B kxÿ ykV .

We will prove the following proposition.

Theorem 3.6. For any m P 1, Sm admits a compact attractor Am which attracts bounded set of Vm. Moreover,dHÿ1�Am;A� ! 0 as m!1.

This proposition is a consequence of the Theorem 3.5 with h � 1=m. Let us prove that S and Sh verify thehypothesis of Theorem 3.5. We take V � Hÿ1�X�; Vh � Vm � Hÿ1�X�; N1 � B�0;R2�, where R2 is a non-negative real number R2 > K2, K2 is de®ned in [9]. Then, using Theorem 3.1, we obtain that for

516 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 11: On the “genetic memory” of chaotic attractor of the barotropic ocean model

N2 � B�0;R2 � K2� the ®rst and the second items of Theorem 3.5 hold for any t0 2 R�. For any d0 > 0 letN3 � B�0;R2 �

���2p

K2 � d0�, then the third item of the Theorem 3.5 is true for T �xh� � �1. Moreover, Sh�t�is asymptotically smooth since the ball B�0;R2�, R2 > K2 attracts any bounded set in Hÿ1�X�.

Let us now prove that Sh�t� conditionally approximates S�t� on N3 uniformly on compact set of R�. Weintroduce the following variational formulation of the problems (33) and (35). Find x�t� 2 Hÿ1�X� suchthat

oxot;u

� �� ��J�w;x� `�;u�� � lhDx;ui ÿ r��x;u�� ÿ ��f ;u�� for 8u 2 Hÿ1�X�;

��x;u�� � ÿ ��rw;ru�� 8u 2 Hÿ1�X�; �36�

where h�; �i is a duality D0�X� � d�X� and ���; ��� is a scalar product de®ned previously in Hÿ1�X�.The problem (35) will be expressed as follows: ®nd xm�t� 2 Vm such as for 8um 2 Vm

dxm

dt;um

� �� ��J�w;x� `�;um�� � l��Dxm;um�� ÿ r��xm;um�� ÿ ��fm;um��;

��xm;um�� � ��rwm;um��: �37�

For any u�t� 2 Hÿ1�X� we denote by �u�t� the projection of u on Vm de®ned by ��uÿ �u;um�� � 0 8um 2 Vm.Since Vm � V , we can write (36) with um 2 Vm. Then, subtracting (37) from (36) we obtain, for any um 2 Vm

o�xÿ xm�ot

;um

� �� ��J�w;x� `�;um�� � ��J�wm;xm � `�;um��

� lhD�xÿ xm�;umi ÿ r���xÿ xm�;um�� ÿ ��f ÿ fm;um��:

We introduce �x 2 Vm such that ��xÿ �x;um�� � 0 and �w 2 Vm such that D �w � �x. Note that, by the equationof Vm, ��wÿ �w;um�� � 0 for any um 2 Vm. So, we obtain for um � �xÿ xm 2 Vm

dk �xÿ xmk2

dt� o�xÿ �x�

ot; �x

ÿ xm

!!� l

2kr �xÿrxmk2 � r

2k �xÿ xmk2

ÿ l��D�xÿ �x�; �xÿ xm�� � ��J�w;xÿ �x�; �xÿ xm�� � ��J�w; �xÿ xm�; �xÿ xm��� ��J�wÿ �w;xm � `�; �xÿ xm�� � ��J� �wÿ wm;xm � `�; �xÿ xm��� ��f ÿ fm; �xÿ xm��: �38�

Now, since o �x=ot � ox=ot we have ���oxÿ �x=ot�; �xÿ xm�� � 0. Moreover, ��D�xÿ �x�; �xÿ xm�� ���xÿ �x;D� �xÿ xm��� � 0 since D� �xÿ xm� 2 Vm. The de®nition of the scalar product in Hÿ1�X� and anintegration by parts show that ��J� �wÿ wm;xm � `�; �xÿ xm�� � 0.

So we have

dk �xÿ xmk2

dt� lkr� �xÿ xm�k2 � rk �xÿ xmk26 j��J�w;xÿ �x�; �xÿ xm��j

� j��J�w; �xÿ xm�; �xÿ xm��j � j��J�wÿ �w;xm � `�; �xÿ xm��j � j��f ÿ fm; �xÿ xm��j: �39�Let us now bound the terms above. Using (3) we get

j��J�wÿ �w; `�; �xÿ xm��j � bowÿ �w

ox; �x

����� ÿ xm

!!����� � bowÿ �w

ox; �w

����� ÿ wm

!�����6 bjr�wÿ �w�jj �wÿ wmj6

b�����k1

p jr�wÿ �w�jjr� �wÿ wm�j

6 rjr� �wÿ wm�j2 �b2����������4rk1

p jr�wÿ �w�j2: �40�

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 517

Page 12: On the “genetic memory” of chaotic attractor of the barotropic ocean model

The Jacobian term can easily be bounded using integration by part and the fact J�f ; g� � 0 on theboundary for all f ; g 2 H 1

0 �X�j��J�w;xÿ �x�; �xÿ xm��j � j�J�w; �wÿ wm�;xÿ �x�j � j�rJ�w; �wÿ wm�;r�wÿ �w�:�j6 krwkH2 jr� �wÿ wm�jL4 jr�wÿ �w�jL4 � jrwjL4k �wÿ wmkH2 jr�wÿ �w�jL4 : �41�

Since w; �w;wm 2 H 10 �X� \ H 2�X�, then kwkH2 6 cjDwjL2 and k �wÿ wmkH2 6 cj �xÿ xmjL2 . Moreover

jf jL4 < cjf j1=2

L2 jrf j1=2

L2 . Then

j��J�w;xÿ �x�; �xÿ xm��j6 cjxjjr� �wÿ wm�j1=2j �xÿ xmj1=2jr�wÿ �w�j1=2jxÿ �xj1=2

� cjrwj1=2jxj1=2jr�wÿ �w�j1=2jxÿ �xj1=2j �xÿ xmj1=2

6 cjxj2jr� �wÿ wm�jjr�wÿ �w�j1=2

� cj �xÿ xmjjxÿ �xjjr�wÿ �w�j1=2

� cjrwjjxjjr�wÿ �w�jjxÿ �xj � lej �xÿ xmj2

6 cjxj2jr� �wÿ wm�j2 � cjxj2jr�wÿ �w�j2

� 2lej �xÿ xmj2 � cjxÿ �xj2jr�wÿ �w�j2

� cjrwjjxjjr�wÿ �w�jjxÿ �xj

6 2lej �xÿ xmj2 � cjxj2jr� �wÿ wm�j2

� c jxj2�

� jxÿ �xj2 � jrwjjxjjxÿ �xj�jr�wÿ �w�j:

Notice that jxÿ �xj2 � ��xÿ �x;D�xÿ �x��� � ��xÿ �x;Dx�� since D �x 2 Vm. Then jxÿ �xj2 �jxj2 ÿ � �x;x�. But �xÿ �x; �x� � 0 so � �x;x� � j �xj2. Then we have jxÿ �xj2 � jxj2 ÿ j �xj26 jxj2. We con-clude that, for any e > 0,

j��J�w;xÿ �x�; �xÿ xm��j6 2lej �xÿ xmj2cjxj2jr� �wÿ wm�j2

� cjxj2 1� � jrwj�jr�wÿ �w�j:The next term in the Eq. (39) can be bounded as

j��J�w; �xÿ xm�; �xÿ xm��j � j�J�w; �xÿ xm�; �wÿ wm�j� j�J�w; �wÿ wm�; �xÿ xm�j6 cj �xÿ xmjrwjL4 jjr� �wÿ wm�jL4

6 cj �xÿ xmjrwj1=2jxj1=2jr� �wÿ wm�j1=2j �xÿ xmj1=2

6 Mej �xÿ xmj2 � cjrwjjxjj �xÿ xmjjr�wÿ wm�j

6 2Mej �xÿ xmj2 � cjrwj2jxj2jr� �wÿ wm�j2:

And the next one

j��J�wÿ �w;xm�; �xÿ xm��j � j�J�wÿ �w; �wÿ wm�;xm�j6 cjxmjjr�wÿ �w�jL4 jr� �wÿ wm�jL4

6 cjxmjjr�wÿ �w�j1=2jxÿ �xj1=2jr� �wÿ wm�j1=2j �xÿ xmj1=2

6 cjxmjjr�wÿ �w�jjxÿ �xj � cjxmjjr� �wÿ wm�j �xÿ xmj6 j �xÿ xmj2 � cjxmj2jr� �wÿ wm�j2 � cjxmjjxjjr�wÿ �w�j:

518 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 13: On the “genetic memory” of chaotic attractor of the barotropic ocean model

Source term in the Eq. (39) will be written as

��f ÿ fm; �xÿ xm��6 kf ÿ fmkÿ1k �xÿ xmkÿ1

6 1

2kf ÿ fmk2

ÿ1 �1

2k �xÿ xmk2

ÿ1 �1

2kf ÿ fmk2

ÿ1 �1

2jr� �wÿ wm�j2:

Inserting all these results in (39) we obtain replacing the Hÿ1 norm by the L2 one

1

2

djr� �wÿ wm�j2dt

� lj �xÿ xmj2 � rjr� �wÿ wm�j2

6 rjr� �wÿ wm�j2 � cjr�wÿ �w�j � 5lej �xÿ xmj2

� cjxj2r� �wÿ wm�j2 � cjxj2 1� � jrwj�jr�wÿ �w�j� cjrwj2jxj2jr� �wÿ wm�j2 � cjxmj2jr� �wÿ wm�j2

� cjxmjjxjjr�wÿ �w�j2

� 1

2kf ÿ fmk2

ÿ1 �1

2jr� �wÿ wm�j2 �42�

then for e � 5

1

2

djr� �wÿ wm�j2dt

6 cjr�wÿ �w�j2 � cjxj jxj� � jxjjrwj � jxmj�jr�wÿ �w�j

� c jxj2�1�

� jrwj� � jxmj2 � 1

2

�� jr� �wÿ wm�j2

� 1

2kf ÿ fmk2

ÿ1: �43�

We recall that jr�wÿ �w�j2 � jrwj2 ÿ j �wj26 jrwj2 and by the Theorem 3.1 jrw�t�j26 jrw�0�j2 � K 6R forw�0� 2 N3. Then

1

2

djr� �wÿ wm�j2dt

6 c jxj2�

� jxmj2 � 1�jr� �wÿ wm�j2

� c jxj2�

� jxjjxmj � 1�jr�wÿ �w�j � 1

2kf ÿ fmk2

ÿ1:

Now, by de®nition of �w we have

jr�wÿ �w�j � kxÿ �xk2ÿ1 � ��xÿ �x;xÿ vm � vm ÿ �x�� 8vm 2 Vm

6 ��xÿ �x;xÿ vm�� since vm ÿ �x 2 Vm

6 kxÿ �xkÿ1kxÿ vmkÿ1 8vm 2 Vm �44�so

jr�wÿ �w�j6 infvm2Vmkxÿ vmkÿ1 �45�

and

djr� �wÿ wm�j2dt

6 c jxj2�

� jxmj2 � 1�jr� �wÿ wm�j2 � kf ÿ fmk2

ÿ1

� c jxj2�

� jxjjxmj � 1�

supt2�0;T �

infvm2Vmkxÿ vmkÿ1: �46�

Integrating in time between 0 and T and using the Gronwall lemma we get

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 519

Page 14: On the “genetic memory” of chaotic attractor of the barotropic ocean model

jr� �w�t� ÿ wm�t��j26 jr� �w�0� ÿ wm�0��j2 � Tkf ÿ fmk2

ÿ1

� c supt2�0;T �

infvm2Vmkxÿ vmkÿ1

Z T

0

�jxj2 � jxjjxmj � 1�dt

!

� exp

Z T

0

�jxj2�

� jxmj2 � 1�dt�: �47�

Finally, using the result of the Theorem 3.1 we conclude that for 8t 2 �0; T �jr� �w�t� ÿ wm�t��j26K�T ; jrw�0�j�

� jr� �w�0�

ÿ wm�0��j2 � kf ÿ fmk2ÿ1 � sup

t2�0;T �inf

vm2Vmkxÿ vmkÿ1

!: �48�

By de®nition of xm�0� we have r �w�0� � rwm�0�, so

jr� �w�t� ÿ wm�t��j26K�T ; jrw�0�j� supt2�0;T �

infvm2Vmkx

ÿ vmkÿ1 � kf ÿ fmkÿ1

!: �49�

Then

jr�wÿ wm�j6 jr�wÿ �w�j � jr� �wÿ wm�j

6K 0�T ; jrw�0�j� supt2�0;T �

infvm2Vmkx

ÿ vmk1=2

ÿ1 � kf ÿ fmkÿ1

!:

And ®nally

jr� �w�t� ÿ wm�t��j26 g�T ; jrw�0�j;m�; �50�where g�T ; jrw�0�j;m� ! 0 when m!1. The proof is completed.

4. Experiment with low dissipation

We perform two experiments for di�erent values of forcing and friction parameters. So far we intend tostudy chaotic behaviour of the ¯ow, we need to choose parameters to have a developed turbulence.

It is known [25], there exist two types of transition to instability in models of this kind. The major typesof destabilisation are identi®ed to be either the meandering of the mid-latitude eastward jet, or the Rossbywave radiation from the westward return ¯ow. In this paper we want to realize both destabilisationmechanisms. In order to do this we perform two experiments with di�erent dissipation parameters. It turnsout to be su�cient to double the friction in order to get another type of destabilisation.

A ®rst approach of the wind driven ocean is provided by the Sverdrup's balance between b e�ect andwind forcing, b�ow=ox� ' f giving a maximum y velocity component vS � �kf k=b� � �2ps0=bqHL�.

In all experiments we use several values of the wind tension from the range de®ned in Table 1. Thisvalues range provides di�erent type of solution of the model from a stable stationary point to chaotic

Table 1

Parameters

Parameter Experiment 1 Experiment 2

Bottom drag r 5� 10ÿ8 sÿ1 10ÿ7 sÿ1

Lateral friction l 1250 m2=s 2500 m2=s

Wind stress range s0 0:6 . . . 1:1 dyne=cm2 1:0 . . . 3:25 dyne=cm2

Corresponding Sverdrup velocity vS 1 . . . 1:7 cm=s 1:6 . . . 5:1 cm=s

520 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 15: On the “genetic memory” of chaotic attractor of the barotropic ocean model

behaviour. The upper limit of this range is determined as the largest value for which the scheme remainsstable.

Our purpose is to identify the characteristic features of the chaotic behaviour of the model and todetermine the origins of these features basing on the hypothesis that they realize the ``memory'' ofchaos.

4.1. Chaotic solution

Largest forcing in the ®rst experiment produces irregular statistically stable solution. The time averagestreamfunction presented in Fig. 2A is antisymmetric with respect to the middle of basin. The variability ofthe streamfunction is illustrated by its dispersion

g2�x; y� � 1

T

Z T

0

�w�x; y; t� ÿ �w�x; y��2 dt where �w�x; y� � 1

T

Z T

0

w�x; y; t�dt: �51�

To distinguish the principal modes of the variability in this system we perform the EOF analysis. We lookfor an orthogonal basis in the phase space of the model in sense of the scalar product (6). We require the®rst vectors contain maximum of dispersion. This basis is can be found as the eigen basis of the covariancematrix.

C/ � k/; Ci;j � 1

T

Z T

0

�Mw�iwj dt: �52�

Eigenvalues kn show the dispersion associated with the eigenvector /n. We arrange the sequence ofeigenvectors, or EOFs so that their corresponding eigenvalues decrease. So the largest dispersion is asso-ciated with the ®rst EOFs.

One can see in Fig. 3A that the principal dispersion is associated with the perfectly symmetric ®rst EOFcomponent.

To distinguish what kind of processes are dominated we construct the probability density function(PDF) of the projection of the streamfunction on the ®rst EOF vector. One can clearly see in Fig. 3B thebimodality of the distribution. That means the system spend much time in a state with high projection on

Fig. 2. (A) Mean streamfunction in the experiment 1 with s0 � 1:1 dyne=cm2 and (B) Streamfunction dispersion in the experiment 1

with s0 � 1:1 dyne=cm2.

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 521

Page 16: On the “genetic memory” of chaotic attractor of the barotropic ocean model

the ®rst EOF, either positive or negative. Hence, the system spend the major time in states which possessesno symmetry.

These bimodality allows us to speak on the existence of two basins of the global attractor, each char-acterised by the sign of the projection on the ®rst EOF vector, either positive or negative. The mean lifetimeof the solution in each basin of attractor is about 300 yr.

To verify the transitions of the solution between basins occur due to the dynamics rather than due tosome numerical noise we add an arbitrary small forcing to the right hand side of the model. This additionalforcing is a ``white noise'' distributed of an amplitude 0:01� kf k. The addition does not change the fre-quency of transition and the mean lifetime in basins of attraction remains the same.

Thus, the system spends the major time in one of two clusters. Transitions between clusters are irregular.They occur once a �300 yr, i.e. we can say the system is almost intransitive.

We reconstruct the cluster centres as a linear combination of the mean state Fig. 2A and the ®rst EOFFig. 3A normalised by its dispersion square root. This state is presented in Fig. 4A.

One can easily see two clear maxima in the Energy spectrum presented in Fig. 4B. These maximacorrespond to oscillations with periods 75 and 150 days respectively.

To verify if it is a chaotic solution and to estimate its attractor dimension we calculate its globalLyapunov exponents (21). To be sure the composition length chose is su�ciently large to approximate thein®nity limit, we calculate averaged local exponents for di�erent composition lengths. The interval of 8960days has been divided into 26ÿL subintervals of 2L � 140 days for each L in 06 L6 6. Local exponents have

Fig. 3. (A) First EOF (k1 � 6:5� 107) in the experiment 1 with s0 � 1:1 dyne=cm2 and (B) PDF of the streamfunction projection on

the ®rst EOF (k1 � 6:5� 107) in the experiment 1 with s0 � 1:1 dyne=cm2.

Fig. 4. (A) Mean state +�����k1

p � EOF in the experiment 1 with s0 � 1:1 dyne=cm2 and (B) Energy spectrum in the experiment 1 with

s0 � 1:1 dyne=cm2.

522 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 17: On the “genetic memory” of chaotic attractor of the barotropic ocean model

been calculated for each subinterval and averaged over all subintervals of this composition length. Thus weget averaged local exponents as a function of their composition length in the range 140 . . . 8960 days, whichcorrespond to the interval 06 L6 6.

In order to increase the accuracy of estimates of global exponents, local values have been extrapolated toa longer time by the relation (28). The value k1 found by application (28) we consider as the global ex-ponent value Fig. 5.

So far there exist two positive Lyapunov exponents we deal with a chaotic behaviour of the solution.Among the principal features of this chaos we can note· long-time integration provides symmetrical average solution as well as its dispersion,· the ®rst EOF is symmetrical with a maximum near the jet,· the distribution function of the projection on the ®rst EOF is bimodal. The solution spent the major time

in one of two clusters with very rare transitions between them.

Fig. 5. Local Lyapunov exponents for di�erent composition length. Experiment 1 with s0 � 1:1 dyne=cm2. Limit values:

k1 � 0:078; k2 � 0:015. Attractor dimension� 5.8.

Fig. 6. (A) Mean streamfunction in the experiment 1 with s0 � 0:6 dyne=cm2 and (B) Eigenmode of the linearised operator in the

experiment 1 with s0 � 0:6 dyne=cm2.

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 523

Page 18: On the “genetic memory” of chaotic attractor of the barotropic ocean model

· There are two maxima in the energy spectrum. They correspond to oscillations with periods 75 and 150days respectively.To ®nd the origins of this chaos and explain its characteristic features we shall look for a possible bi-

furcations which may occur during the period of transition to chaos. We choose in particular the Grasgofnumber which is proportional to the amplitude of the forcing as the bifurcation parameter.

4.2. Symmetric stationary solution and Pitchfork bifurcation

When forcing is lower than 0:6 dyne=cm2 there exists a symmetric stable stationary solution Fig. 6A.The eigenmode of the operator linearised with respect to this solution with largest real part of eigenvalueis real Fig. 6B. When it looses its stability, a couple of asymmetric stationary solutions arises in thePitchfork bifurcation. According to the Theorem 3.2 the total number of stationary points must be odd,

Fig. 7. (A) Mean streamfunction in the experiment 1 with s0 � 0:8 dyne=cm2 and (B) Module of the eigenmode of the linearised

operator in the experiment 1 with s0 � 0:8 dyne=cm2.

Fig. 8. (A) Mean streamfunction in the experiment 1 with s0 � 0:85 dyne=cm2 and (B) Streamfunction dispersion in the experiment 1

with s0 � 0:85 dyne=cm2.

524 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 19: On the “genetic memory” of chaotic attractor of the barotropic ocean model

we can conclude the symmetric stationary point does not disappear. It remains stationary while becomesunstable.

The global attractor of the system in this case is a connected set composed of the three stationary points.Instable manifold of the symmetric point is included in the attractor [24]. So far the only real eigenvalue ofthe linearised operator becomes positive in the vicinity of the bifurcation point, we can state that the di-mension of the global attractor is equal to 1. The attention should be focussed on the fact that the di-mension of the global attractor is not equal to the dimension of the union of local attractors, which is theset of two stable equilibria in this case.

So we can conclude that after the bifurcation the system becomes intransitive with two invariant localattractors represented by two points.

The eigenmode of the tangent operator linearised with respect to this equilibrium is a real mode. Thismode is shown in Fig. 6B.

Fig. 9. (A) Energy spectrum. Experiment 1 with s0 � 0:85 dyne=cm2 and (B) Energy spectrum. Experiment 1 with

s0 � 0:87 dyne=cm2.

Fig. 10. (A) Mean streamfunction in the experiment 2 with s0 � 3:25 dyne=cm2 and (B) Streamfunction dispersion in the experiment 2

with s0 � 3:25 dyne=cm2.

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 525

Page 20: On the “genetic memory” of chaotic attractor of the barotropic ocean model

One can easily see that the equilibrium of the model with low forcing is perfectly symmetrical as well asthe average state of the long-time integration experiment Fig. 2A. And the ®rst eigenmode of the tangentlinear operator is similar to the ®rst EOF vector obtained in the long-time integration Fig. 3A.

4.3. Asymmetric stationary solution and Hopf bifurcation

Two asymmetric stationary solution exist in the range of the wind stresses from 0.6 to 0.85 dyne=cm2.Streamfunction pattern of one of them is presented in Fig. 7A. The second one can be obtained by sym-metric transformation with respect to the middle latitude.

The square of the eigenmode of the tangent operator linearised with respect to this equilibrium is shownin Fig. 7B.

At the wind tension 0:85 dyne=cm2 we observe the Hopf bifurcation. Each stationary solution Fig. 7Alooses its stability and produce a periodic orbit. The eigenmode of the linearised operator has complexconjugate pair of eigenvalues with imaginary part 0.087 dayÿ1 that corresponds to a period of about 75days The square of this eigenmode is shown in Fig. 7B.

The local attractor dimension after the bifurcation becomes equal to 1. The global attractor, beingcomposed of (at least) three stationary points and two limited loops has a dimension greater or equal to2. And again in this case, the system is also intransitive, because the global attractor includes two localones.

4.4. Periodic solution and Feigenbaum phenomenon

When the wind tension is comprised between 0.85 and 0.90 dyne=cm2 we observe the Feigenbaumphenomenon of doubling period. At the forcing value 0.85 the solution is stable periodic with period2p=0:087 � 72 days. Mean streamfunction and its dispersion are shown in Fig. 8.

One can ®nd that the dispersion of this periodic solution is identical to the eigen mode of the linearisedoperator of the model with the wind tension s0 � 0:8 dyne=cm2. This fact, together with the equality be-tween the period of the periodic solution and the imaginary part eigen mode of the linearised operator,con®rms the Hopf bifurcation really took place.

At the forcing value 0.88 one can observe the ampli®cation of the double period and quadruple periodmodes. This fact is illustrated by the energy spectrum of the solution in Fig. 9.

Hence, the chaotic attractor of the model has been formed by the sequence of two symmetric doublingperiod bifurcation sets resulted in two symmetric basins of global attractor. These basins are linked by a``very weak'' connection providing very rare transitions between basins.

5. Experiment with high dissipation

Analysis of unstable eigenmodes of the linearised operator shows that transformations of energy fromthe mean ¯ow to the oscillations occurs near the jet, in the region of high velocity gradients. However, as ithas been remarked above, there exists another possible mechanism of the transformation: emission of theRossby waves from the Eastern boundary. To realize this mechanism we increase the dissipation and theforcing in the model.

So, the second experiment was performed for the model with the lateral friction and the bottom drag twotimes higher.

5.1. Chaotic solution

Largest forcing in Table 1 in the second experiment produces irregular statistically stable solution. Thetime average streamfunction presented in Fig. 10A is antisymmetric with respect to the middle of basin. Thevariability of the streamfunction is illustrated by its dispersion (51).

526 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 21: On the “genetic memory” of chaotic attractor of the barotropic ocean model

Along with the maximum near the jet, one can see two clear maxima near the western boundary in Fig.10B. This con®rms that the Rossby wave emission plays an important role in this experiment.

To distinguish the principal modes of the variability in this system we perform the EOF analysis as wellas in the experiment 1.

One can easily see that, similarly to the experiment 1, the principal dispersion is associated with theperfectly symmetric ®rst EOF component. The PDF function of the projection of the streamfunction on the®rst EOF vector has also two maxima as in the experiment 1. One can clearly see the bimodality of thisdistribution in Fig. 11B. That means the system spends much time in a state with high projection on the ®rstEOF, either positive or negative.

To verify it is a chaotic solution with a strange attractor we calculate its Lyapunov exponents (21) andwe perform an extrapolation of mean local exponents to in®nite time by the relation (28).

So far there exist four positive Lyapunov exponents we deal with a chaotic behaviour of the solution. To®nd the origins of this chaos we shall look for a possible bifurcations which may occur during the transitionperiod to chaos. We choose the amplitude of the forcing as the bifurcation parameter.

One can easily see one clear maximum in the energy spectrum presented in Fig. 12B. This maximumcorresponds to the Rossby waves with period 17 days.

Among the principal features of the chaotic attractor in this case we can note· long-time integration provides symmetrical average solution as well as its dispersion,· the ®rst EOF is symmetrical with a maximum near the jet,· the distribution function of the projection on the ®rst EOF is bimodal. The solution spent the major time

in one of two clusters with rare transitions between them.

Fig. 11. (A) First EOF (k1 � 8:8� 108) in the experiment 2 with s0 � 3:25 dyne=cm2 and (B) Distribution of the projection of the

streamfunction on the ®rst EOF (k1 � 8:8� 108) in the experiment 2 with s0 � 3:25 dyne=cm2.

Fig. 12. (A) Local Lyapunov exponents for di�erent composition length. Experiment 2 with s0 � 3:25 dyne=cm2. Limit values:

k1 � 0:095; k2 � 0:057. Attractor dimension� 9.8 and (B) Energy spectrum in the experiment 2 with s0 � 3:25 dyne=cm2.

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 527

Page 22: On the “genetic memory” of chaotic attractor of the barotropic ocean model

· There is one maximum in the energy spectrum. This maximum corresponds to oscillations with period 17days.

5.2. Symmetric stationary solution and Pitchfork bifurcation

When forcing is lower than 1:2 dyne=cm2 there exists a symmetric stable stationary solution. Theeigenmode of the operator linearised with respect to this solution with largest real part of eigenvalue is realFig. 13B.

One can easily see in Fig. 13A that the equilibrium of the model with low forcing is perfectly symmetricalas well as the average state of the long-time integration experiment Fig. 10A. And the ®rst eigenmode of thetangent linear operator in Fig. 13B is similar to the ®rst EOF vector obtained in the long-time integrationFig. 11A.

Fig. 13. (A) Mean streamfunction in the experiment 2 with s0 � 1:0 dyne=cm2 and (B) Eigenmode of the linearised operator in the

experiment 2 with s0 � 1:0 dyne=cm2.

Fig. 14. (A) Stationary streamfunction in the experiment 2 with s0 � 2:1 dyne=cm2 and (B) Module of the eigenmode of the linearised

operator in the experiment 2 with s0 � 2:1 dyne=cm2.

528 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 23: On the “genetic memory” of chaotic attractor of the barotropic ocean model

When the wind tension reaches the value 2:1 dyne=cm2 the solution Fig. 13A looses its stability, a coupleof asymmetric stationary solutions arises in the Pitchfork bifurcation. As well as in the experiment 1, ac-cording to the theorem Theorem 3.2 the total number of stationary points is odd, the symmetric stationarypoint does not disappear. It remains stationary while becomes unstable. The global attractor of the systemin this case is a connected set composed of the three stationary points.

So we can conclude as in the experiment 1, that after the bifurcation the system becomes intransitive withtwo invariant local attractors represented by two points.

5.3. Asymmetric stationary solution, direct and inverse Hopf bifurcations

Stationary solution Fig. 14A looses its stability and generate a periodic orbit. The eigenmode of thelinearised operator has complex conjugate pair of eigenvalues with imaginary part 0.244 dayÿ1 that cor-responds to a period of about 25:7 days. The square of this eigenmode is shown in Fig. 14B.

Fig. 15. (A) Mean streamfunction in the experiment 2 with s0 � 2:2 dyne=cm2 and (B) Streamfunction dispersion in the experiment 2

with s0 � 2:2 dyne=cm2.

Fig. 16. (A) Stationary streamfunction in the experiment 2 with s0 � 2:4 dyne=cm2 and (B) Module of the eigenmode of the linearised

operator in the experiment 2 with s0 � 2:4 dyne=cm2.

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 529

Page 24: On the “genetic memory” of chaotic attractor of the barotropic ocean model

At the wind tension is comprised between 2.1 and 2.8 we observe the sequence of three Hopf bifurca-tions. The ®rst one occurs at the forcing 2.1 dyne=cm2.

At the forcing value 2.2 the solution is stable periodic with period 2p=0:244 � 25:7 days. Meanstreamfunction and its dispersion are shown in Fig. 15.

One can ®nd that the dispersion of this periodic solution is identical to the eigenmode of the linearisedoperator of the model with the wind tension s0 � 1:0 dyne=cm2. This fact, together with the equality be-tween the period of the periodic solution and the imaginary part eigenmode of the linearised operator,con®rms the Hopf bifurcation really took place.

However, instead the Feigenbaum phenomenon as in experiment 1, further increasing of the forcingvalue results in the inverse Hopf bifurcation. The periodic solution become stationary and stable. Thestreamfunction pattern and the module of eigenmode for the model with forcing 2:4 dyne=cm2 are shown inFig. 16.

5.4. Periodic solution and transition to chaos

When the wind tension is comprised between 2.8 and 3.0 dyne=cm2 we observe the transition to chaoswhich is not realized by the Feigenbaum phenomenon as before. While at forcing 2.8 the solution is stableperiodic with period 2p=0:380 � 16:5 days. At forcing 2.85 one can see ampli®cation of harmonics withperiods 2p=0:305 � 20:5 days and 2p=0:060 � 104 days. Comparing their frequencies with imaginary partsof eigen modes of the linearised operator, one can conclude that the reason of this ampli®cation is thesecond and third Hopf bifurcation.

One can see energy spectra of the solution with forcing 2:8 and 2:85 dyne=cm2in Fig. 17.Hence the chaotic attractor of the model has been formed by the sequence of two symmetric bifurcation

sets resulted in two symmetric basins of global attractor. These basins are linked by a ``weak'' connectionproviding very rare transitions between basins.

Fig. 17. (A) Energy spectrum. Experiment 2 with s0 � 2:80 dyne=cm2 and (B) Energy spectrum. Experiment 2 with

s0 � 2:85 dyne=cm2.

530 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532

Page 25: On the “genetic memory” of chaotic attractor of the barotropic ocean model

The transition to chaos has not been realized by Feigenbaum phenomenon of doubling period, so we®nd the only clear maximum in the energy spectrum of the chaotic solution Fig. 12B which corresponds tothe principal frequency of the periodic solution Fig. 17.

6. Conclusion

It is shown in this paper that the structure of attractor of barotropic ocean model can be partially ex-plained by the sequence of bifurcations that the system is subjected by variations of the leading parameters.The principal feature of the studied system seems to be the rise of the couple of non-symmetric stablestationary solutions in the model with symmetric forcing. This couple arises in the bifurcation of thesymmetric stable equilibrium. The global attractor of the system includes two local attractors and thereforethe system is intransitive. This fact provides the existence of two ``almost invariant'' basins of chaoticattractor with very rare transitions between them (average lifetime in a basin is about 300 yr in the ®rstexperiment and about 10 yr in the second one). The existence of these basins is illustrated by the bimodalityof the PDF of the ®rst EOF component of the chaotic solition.

The ``memory'' of chaos appears also in the presence of maxima in the spectrum of energy. Thesemaxima correspond either to the principal frequency of the limit cycle which arose in the Hopf bifurcation,or to the frequencies of the Feigenbaum phenomenon. The physical reason of the principal oscillations canbe di�erent, either emission of the Rossby waves by the Eastern boundary, or nonlinear interaction of thejet with the Rossby waves.

The results of this work would be more reliable to the real processes in the ocean when the morecomprehensive ocean model with a high spatial resolution is used. This work is planned to the nearestfuture.

Acknowledgements

This work has been performed in frames of project of the French±Russian Lyapunov Institute. Theprincipal part of this work has been carried out during the visit of V. Dymnikov in the Institute of ElieCartan, University of Nancy 1 thanks to the support of the French Ministry of High Education and Re-search. All the contour pictures have been prepared by the Grid Analysis and Display System (GrADS)developed in the Centre for Ocean±Land±Atmosphere Interactions, Department of Meteorology, Uni-versity of Maryland.

References

[1] V. Dymnikov, A. Filatov, Stability of the Large-scale Atmospheric Processes, Hydrometeoizdat, Leningrad, 1990.

[2] V. Dymnikov, A. Filatov, Mathematics of Climate Modelling, Birkhauser, Basel, 1996.

[3] V.P. Dymnikov, E.V. Kazantsev, On the attractor structure generated by the system of equation of the barotropic atmosphere,

Izvestiya, Atmospheric and Oceanic Physics 29 (5) (1993).

[4] V. Dymnikov, A. Gritsun, On the structure of the attractors of the ®nite-dimensional approximations of the barotropic vorticity

equation on the rotating sphere, Russ. J. Numer. Anal. Math. Modelling 12 (1) (1997) 13±32.

[5] Ch. Bernier, E. Kazantsev, Theoretical and numerical estimates of instantaneous Lyapunov exponents of quasi-geostrophic ocean

dynamics, Appl. Math. Comp. Sci. 6 (2) (1996) 223±245.

[6] S. Jiang, F.F. Jin, M. Ghil, Multiple equilibria, periodic and aperiodic solutions in a wind-driven double gyre shallow-water

model, J. Phys. Oceanography 25 (1995) 764±785.

[7] H. Itih, M. Kimoto, Multiple attractors and chaotic itinerancy in a quasi-geostrophic model with realistic topography:

implications for weather regimes and low-frequency variability, J. Atm. Sci. 53 (15) (1996) 2217±2231.

[8] B. Wang, Z. Fang, Chaotic oscillations of tropical climate: a dynamic system theory for ENSO, J. Atm. Sci. 53 (19) (1996) 2786±

2802.

[9] Ch. Bernier, Existence of attractor for the quasi-geostrophic approximation of the Navier±Stockes equations and estimates of its

dimension, Adv. Math. Sci. Appl. 4 (2) (1994) 465±489.

[10] A.A. Il'in, The Navier±Stokes and Euler equations on two-dimensional closed manifolds, Math. USSR Sbornik 69 (1) (1991).

V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532 531

Page 26: On the “genetic memory” of chaotic attractor of the barotropic ocean model

[11] A.N. Filatov, Estimation of the number of instable stationary solutions for spectral models of the atmosphere, Works of the

Russian Hydro-meteorological Centre 323 (1992).

[12] W.R. Holland, The role of mesoscale eddies in the general circulation of the ocean-numerical experiments using wind driven quasi-

geostrophic model, J. Phys. Oceanography 8 (3) (1978).

[13] G.J. Fix, Finite elements models for ocean circulation problems, J. Appl. Math. 29 (1975).

[14] C. Le Provost, C. Bernier, E. Blayo, A comparison of two numerical methods for integrating a quasi-geostrophic multilayer model

of ocean circulations: ®nite element and ®nite di�erence methods, J. comput. phys. 110 (2) (1994).

[15] M. Bernadou, Modulef: une biblioth�eque modulaire d'�el�ements ®nis, INRIA (1988).

[16] Ch. Bernier, Etude et parall�elisation d'un code d'�el�ements ®nis pour la mod�elisation quasi-g�eostrophique des circulations

oc�eaniques, Ph.D. thesis of INPG, Grenoble, France, 1990.

[17] V.I. Oseledets, A Multiplicative Ergodic theorem. Lyapunov characteristic numbers for Dynamical Systems, Moscow Math. Soc.

19 (1968) 197.

[18] J.L. Kaplan, J.A. Yorke, Chaotic behaviour in multidimensional di�erence equations, Lecture notes in mathematics 730 (1979).

[19] H.D.I. Abarbanel, R. Brown, M.B. Kennel, Lyapunov exponents in chaotic systems: their importance and their evaluation using

observed data, Invited Rev. Modern phys. Lett. B (1991).

[20] E. Kazantsev, Local Lyapunov exponents of the quasi-geostrophic ocean dynamics, Rapport de recherche INRIA 3074 (1996).

[21] V.P. Dymnikov, E.V. Kazantsev, V.V. Kharin, Information entropy and local Lyapunov exponents of barotropic atmospheric

circulation, Izvestiya, Atmospheric and Oceanic Phys. 28 (6) (1992).

[22] S.J. Jacobs, A note on multiple ¯ow equilibria, Pageoph 130 (4) (1989).

[23] Hale, Lin, Raugel, Upper semi-continuity of attractors for approximations of semi-groups and partial di�erential equations,

Math. Comput. 50 (181) (1988) 89±123.

[24] J.P. Eckmann, D. Ruelle, Ergodic theory of chaos and strange attractors, Rev. Modern Phys. 57 (3) (1985) 617±656.

[25] C. Le Provost, J. Verron, Wind-driven ocean circulation transition to barotropic instability, Dyn. Atmos. Oceans 11 (1987).

[26] K. Mo, M. Ghil, Statistics and dynamics of persistent anomalies, J. Atmos. Sci. 44 (1987) 877±901.

[27] K. Mo, M. Ghil, Cluster analysis of multiple planetary ¯ow regimes, J. Geophys. Res. 83 (1988) 10927±10952.

532 V. Dymnikov et al. / Chaos, Solitons and Fractals 11 (2000) 507±532