5 Global dynamical properties

  
5.1 Parameter space and stability of main periodic orbits

Let us now investigate the global dynamical properties of the logarithmic potential in the parameter space (h,q). In what follows we use the formulation given in Sect. 3, denoting the dimensionless energy and time again by h and t, and the variables as lowercase. To get the "true'' energy (i.e., in the old variables) just compute $p_0^2(h+\ln r_{\rm c})$ with p₀²=2 and $r_{\rm c}=0.1$ (see Eq. (17) and discussion around).

First we consider the parameter space. In Sect. 2.2 we show that q should be restricted to some domain, $q_0(h)\le q\le 1$, to have the logarithmic model physical sense; $q^2_0(h)=0.5-\exp(-2h)$. This subspace defines in the hq-plane the "physical region'' where the density is positive everywhere. This is sketched in Fig. 13 where we see that for , the physical region is restricted to . For $h\le \ln 2/2$no bounds for q exist. Therefore to explore the dynamics for generic values of the energy ( ), q₀^*=0.71 is a good lower bound for q.

Figure 13: Parameter space (h,q) of the logarithmic Hamiltonian (17). The physical region is determined by $q_0(h)\le q\le 1$ where $q_0^2(h)=0.5- {\rm e}^{-2h}$ (see Sect. 2.2 and text)

Next we discuss the stability of the long-axis and short-axis periodic orbits. The stability analysis is done by means of a standard Floquet technique and using the fact that both orbits are periodic orbits of 1D Hamiltonians with period T_x and T_y, respectively. Let us summarise the procedure.

Figure 14: a) Stability of the long and short-axis periodic orbits. Horizontal axis: dimensionless energy, h; vertical axis: $\beta$. For the short-axis orbit, $\beta =q^2<1$, while $\beta =1/q^2>1$ corresponds to the long-axis orbit. Grey: instability regions, ${\rm Tr}(M)>2$; white: stability regions, ${\rm Tr}(M)<2$, lines: $\vert{\rm Tr}(M)\vert=2$. b) The same but restricted to the physical domain defined by the curves $q_0^2(h)\le \beta \le 1/q_0^2(h)$. Dark grey: instability "r'' egions inside the physical domain. See text for details

From (17) and denoting, as before, by $\phi(x,y)$ the potential, for the long-axis periodic orbit, we have

$$% \dot{p}_x=-\phi_x(x,0)=-{x\over x^2+1}, \quad \dot{x}=p_x.$$ (27)

Let $\varphi_x(t)$ be a T_x-periodic solution for x(t). Consider now the normal variational equations to the latter vector field for $\varphi_x(t)$, that is,

$$% \dot{\delta p}_y=-\phi_{yy}(\varphi_x(t),0)\delta y=-{1/q^2... ...r 1+ \varphi_x^2(t)}\delta y,\quad \dot{\delta y}=\delta p_y.$$ (28)

In a similar way, for the short-axis periodic orbit, we obtain

$$\begin{eqnarray*}\dot{p}_y=-{y/q^2\over y^2/q^2+1}, \qquad \dot{y}=p_y, \end{eqnarray*}$$

which, after scaling to z=y/q, p_z=qp_y and s=t/q, can be written as

$$% p'_z=-q{z\over z^2+1},\quad z'={p_z\over q},$$ (29)

where prime means derivative respect to s. Denote with $\varphi_z(s)$ a S_y-periodic solution for z(s) (S_y=T_y/q). The normal variational equations are then

$$% \delta p_x'=-{q\over 1+\varphi_z^2(s)} \delta x,\quad \delta x'=q\delta p_x.$$ (30)

The equations of motion (27), (29), and the variational ones (28), (30), can be reduced to a single pair of coupled differential equations of the form

$$% {{\rm d}^2w\over{\rm d}\tau^2}+{w\over 1+w^2}=0,\quad {{\rm d}^2\xi\over{\rm d}\tau^2}+{\beta\over 1+w^2}\xi=0,$$ (31)

where, for the x-axis orbit, $\tau=t$ and $\beta=1/q^2$ while, for the y-axis orbit, $\tau=s=t/q$ and $\beta=q^2.$ Therefore, the stability of both orbits can be studied in the $h\beta$-plane. The region defined by $q_0^2< \beta<1$ corresponds to the y-axis orbit while that where $1<\beta<1/q_0^2$to the x-axis orbit.

The first in (31) leads to

$$% {1\over 2}\left({{\rm d}w\over{\rm d}\tau}\right)^2+{1\over 2} \ln(1+w^2)=h.$$ (32)

From the latter, it is possible to derive asymptotic expressions for the period, P(h), of any orbit at a given h level. Indeed, after change of variables, it is easy to show that

$$\begin{eqnarray*}P(h)=4 {\rm e}^h I(h),\quad I(h)=\int_0^{\sqrt{2h}}{{\rm e}^{-x^2} {\rm d}x\over \sqrt{{\rm e}^{-x^2}-{\rm e}^{-2h}}}, \end{eqnarray*}$$

which for small and large energies reduces to

$$\begin{eqnarray*}P(h)&\approx& 2\pi(1-h)\,{\rm e}^h\quad h\ll 1,\nonumber\\ P(h)&\approx& 2\sqrt{2\pi}\,{\rm e}^h\quad h\gg 1.\nonumber \end{eqnarray*}$$

Here $P=T_x,\, T_y/\sqrt{\beta}$ for the x- and y-axis orbit, respectively. So, in terms of the original variables like in (16), $h\to h/p_0^2-\ln r_{\rm c},\, t\to p_0t/r_{\rm c}$ and for $h\gg p_0^2\ln r_{\rm c}$ we get

$$% T_x(h)\approx{2\sqrt{2\pi}\over p_0}\,{\rm e}^{h/p_0^2},\quad T_y(h)\approx qT_x(h).$$ (33)

In a similar way we have T_x(h) and T_y(h) for small h.

The stability index is defined, as usual, by the trace of the monodromy matrix associated to (31). That is, defining $\eta={\rm d}\xi/{\rm d} \tau$ and $\vec{\xi}=(\xi,\eta)$, the second in (31) can be written as ${\rm d}\vec{\xi}/{\rm d}\tau=D(w;\beta)\vec{\xi}$ where D is the $2\times 2$ differential matrix $D_{11}=D_{22}=0,\, D_{12}=1,\, D_{21}=-\beta/ (1+w^2)$. Thus if we write $\vec{\xi}(\tau)=A(\tau)\vec{\xi}_0$, being $\vec{\xi}_0$ any initial condition and $A(0)={\rm Id}$, then $A(\tau)$satisfies the equation,

$$% {{\rm d}A\over{\rm d}\tau}=D(\varphi_w(\tau);\beta)A$$ (34)

where $\varphi_w(\tau)$ is a P(h)-periodic solution of the first in (31) or of (32). The monodromy matrix is, by definition, M=A(P(h)). To compute M we have to solve (34). However, as the flow is invariant by the change $(\xi,\eta,t)\to (\xi,-\eta,-t)$, then one has A₁₁=A₂₂, using the symplectic character of $A(\tau)$. So to get Tr(M) it is enough to solve for the first column of A. A periodic orbit is linearly stable (unstable) if $\vert{\rm Tr}(M)\vert<2$ (>2). In case $\vert{\rm Tr}(M)\vert=2$, the orbit is said to be marginal stable. For these particular periodic orbits of the logarithmic potential, one has immediately that ${\rm Tr}(M)\ge -2$. So, instability comes only from ${\rm Tr}(M)>2$. This is due to the fact that $D(\varphi_w(\tau))$ is P/2-periodic, so M=A(P/2)² and as A(P/2) is real symplectic, M cannot have negative eigenvalues.

In Fig. 14a we present the stability diagram of both periodic orbits for $0\le h\le 15$ (i.e., a true energy interval about [-4.6,25.4]) and $0\le\beta\le 20$. The region defined by $\beta =q^2<1$ corresponds to the stability diagram of the y-axis orbit while that for $\beta =1/q^2>1$to the x-axis orbit. We plot in grey those zones where ${\rm Tr}(M)>2$, in white those where ${\rm Tr}(M)<2$ and the lines are level curves of $\vert{\rm Tr}(M)\vert=2$. From $(h,\beta)=(0,n^2)$ ($n\ne 0$ integer) emanate the so-called instability tongues. This behaviour is due to the parametric resonance. Indeed, the second in (31) for small h can be written as $\ddot{\xi}+[\beta+\mu(h)c(t)]\xi=0$ where c(t) is periodic and $\mu(0)=0$. For general P-periodic c(t), for instance, $c(t)=\cos(2\pi t/P)$, the instability tongues would emanate from $(\mu,\beta)=(0,n^2/4)$. This is the well known result for the classical Mathieu equation (see, for example, Broer & Levi 1995 and references therein as well as Broer & Simó 1998 for the related unfoldings). However, in our case, the instability tongues emanating from (0,(2k+1)²/4) do not show up and only appear at (0,n²).

It is also worth to mention that the boundaries of the instability zones ( ${\rm Tr}(M)=2$ in the present case) correspond to values of $\beta$ for which the second equation in (31) has a periodic solution with the same period of the first in (31). This periodic solution performs n full oscillations (i.e., it has exactly 2n changes of sign) for any point on the boundary of the n-instability tongue, counted from the bottom. The differences between the two boundaries are as follows: assume that the first in (31) is started, at $\tau=0$ with w=0 and ${\rm d}w/{\rm d}\tau>0$. Then, for the lower boundary, the periodic solution of the variational equation is the first column of the matrix A, the upper boundary is the second column. In an equivalent way, the monodromy matrix, M, has in both cases 1 at the diagonal, and for the lower boundary, M₂₁=0, $M_{12}\ne 0$, while $M_{21}\ne 0$, M₁₂= 0 for the upper one.

Figure 15: Surface of section (x,px), for h=15 and 7 (about 9.8 and 5.8 in dimensionless units) with q=0.95 and 0.75 and $r_{\rm c}=0.1$, p02=2, corresponding to the principal resonances ( $\omega _y/\omega _x=$ a:b) and their associated stochastic layers (see text)

As long as h increases, the unstable zones become wider. Note that the width of these zones approaches asymptotically to a fixed value, $\Delta\beta_n=1$(see Appendix for a proof of this claim). This is more evident for , where the stability regions are extremely narrow at high energies. The convergence of $\Delta\beta_n\to 1$ is the faster the smaller is n. For larger values of $\beta$, , the stability zones (for the x-axis orbit) become significant. Nevertheless almost all of them lie in the unphysical domain. In Fig. 14b we restrict the parameter space to the physical domain, $q_0^2(h)\le \beta \le 1/q_0^2(h)$, where we see that for large h( ), the region is bounded by while for small h ( $h\le\ln 2/2\approx 0.35$), $0<\beta<\infty$. In this figure we plot in dark grey the unstable zones inside the physical region. We observe that the y-axis periodic orbit ($\beta<1$) is always unstable except in some zone of the "harmonic oscillator regime'' ($h\ll 1$). In this regime and for the y-axis orbit is unstable if $q~(=\sqrt{\beta})$ is close to 1, due to the parametric resonance. On the other hand, the x-axis orbit ($\beta>1$) appears to be always stable in the low energy regime. Only a very narrow zone of instability is present, for and ( ). In the range , the x-axis orbit can be either stable or unstable. The stability zone becomes very narrow when we approach $h\approx 7$; $\beta$ is confined to a small interval of the form $(1,1+ \varepsilon(h))$, where $\varepsilon\to 0$ as $\delta\exp(-\gamma h)$, $\gamma$close to 1 and, for instance, $\varepsilon\approx 0.089$, $23.5\ 10^{-6}$, $23.8\ 10^{-10}$ for h=5, 15, 25, respectively. For large energies, , the x-axis orbit is unstable except for q extremely close to 1. Note that for , this orbit is stable in two different $\beta$ intervals, one is that mentioned above and, the second one, starts just beyond the physical region, at $\beta\approx 2$ ( $q\approx\sqrt{2}/2$). The latter result confirms the conjecture given by MS89 about the border of this second stability domain.

To complete this section let us say a few words about the 1:1 periodic orbit. As we have already mentioned in Sect. 2, this orbit is always stable for any physical value of q. Using the formulation given in that section, the rotation period for this orbit can be approximated by

$$% T_{1:1}(h)={2\pi\over\omega_{\theta}}\approx{2\pi\over p_0}\exp \left(-(1+\alpha/2)/2\right){\rm e}^{h/p_0^2},$$ (35)

where $\alpha=(1-q^2)/q^2$. This estimation is obtained using the fact that if $\alpha$ is not too large, the 1:1 orbit is nearly circular and then $T_{1:1}\approx 2\pi a/v_{\rm c}$, where a is the circular radius given by (12) and $v_{\rm c}$, the circular speed, is approximately p₀. From (33) and (35) for $h\gg p_0^2\ln r_{\rm c}$, we obtain

$$\begin{eqnarray*}{T_{1:1}\over T_x}\approx \sqrt{\pi\over 2{\rm e}}{\rm e}^ {-\alpha/4}\approx 0.76\ {\rm e}^{-\alpha/4}, \end{eqnarray*}$$

independent of h. This estimation of the periods agrees, roughly, with the accurate values given by PL96 for relatively small h (-0.4059).

Figure 16: (Top) Width of the stochastic layers where the short and long-axis orbits lie against the dimensionless energy, h, and for three representative values of q. (Bottom)  ${\rm Log}\,\sigma$ after $T(h)\approx 10^4\ T_x(h)$for both layers. See text for details

  
5.2 Irregular motion

In this section we investigate some aspects of the stochastic component of the phase space associated to the logarithmic potential. In Sect. 4.2 we showed, for three values of q (0.9, 0.8, 0.7) and for only one of the energy (-0.4059), that the motion in the logarithmic potential is mostly regular but populated by many resonances. The stochastic component appears to be confined to the gross stochastic layer where the short-axis orbit lies and, in some cases, to the thin layer corresponding to the long-axis orbit. Within the box domain some other chaotic regions are present, for instance, around the 4-periodic island chain that appears close to the gross stochastic layer (see Fig. 9 for q=0.7 and below). To sketch this, in Fig. 15 we show surfaces of section (x,p_x) for comparatively large values of the energy, h=15, 7 and two extreme cases q=0.95, 0.75. In any of these figures we recognize the main resonances, given by $\omega_y/\omega_x$, that lead to loop (1:1), pretzel (4:3), fish (3:2), banana (2:1) and other high order resonances (5:4, 7:5, 5:3, for instance) - see MS89 for the terminology on the periodic orbits and several illustrations. All of them are surrounded by a stochastic layer where the associated unstable periodic orbits lie. It seems that, for these particular values of h and q, all the layers are connected in a single stochastic channel. For comparatively large energies, box orbits disappear leading to boxlets. Certainly, within the loop domain, some isolated chaotic zones should exist.

An interesting aspect to investigate is the connection between different stochastic regions. The results given in Sect. 4.2 for $h\approx -0.4$and q=0.7, show that the stochastic layer around the 1:1 resonance is not connected with that of the 4:3. The same happens with the layers surrounding the 4:3 and 3:2, 3:2 and 5:3, 5:3 and 2:1 resonances. It it clear that the overlap would occur at higher energies. A formal way to investigate connections between different resonances (heteroclinic connections) would be to follow the evolution of any initial condition taken along the unstable manifold associated to each resonance. However, in this case this approach is not so simple. The periodic points associated to, for instance, the 4:3 and the 3:2 resonances change their location and stability properties as h varies. For example, in Fig. 15 the most relevant hyperbolic points lie on the p_x axis, but for lower energies, some of them appear as elliptic on this axis. So, to locate the unstable manifolds, we should look for the hyperbolic points somewhere out of the p_x axis. While this is not a serious complication we have taken a simpler approach.

Therefore we consider the connection between the stochastic layers surrounding the 1:1 and 2:1 resonances. That is, we deal with the problem of finding h such that a point close to the origin at t=0 is close to the border of the section at t=T. So, the derived critical energy, $h_{\rm c}(q),$ is an upper bound to connections with other stochastic layers between them. We study the energy range [-0.4,7.6], about [2.1,6.1] in dimensionless units, with step 0.1 and three representative values of q, 0.9, 0.8 and 0.71. For each energy level we take two orbits such that, at t=0 one is close to (x,p_x)=(0,0) and the other close to (0,p^*_x(h)), where $p_x^*(h)=(2h-p_0^2\ln r_{\rm c}^2)^{1/2}$. Typically we take $(0,\epsilon)$ and $(0,(1-\epsilon)p^*_x(h))$ with $\epsilon=10^{-7}.$ The reason to take $\epsilon$ so small is due to the structure of the stochastic layer (if it exists, that is, if the orbit is unstable). The external part of the layer, close to the border, is populated by many small resonances. If the initial condition falls inside any of them, the orbit will be confined there forever. On the other hand the central part of the layer, around the separatrix, looks like ergodic and therefore it appears to be appropriate to take the initial condition in this part (for details, see Chirikov 1979, Sect. 6).

Next, using the transformation given in Sect. 3, we pass to the variables (x₁, x₃) defined in (20). In this way $0\le \vert x_3\vert\le 1$. Each point is followed during $3.5\ 10^4$ consequents on the Poincaré section y=x₂=0. During the computation we restrict the attention to a small x₁ interval centered at the origin, typically, $\vert x_1\vert\le 0.05$. Within this interval we look for the maximum distance between $(x_1(t),x_3(t))\approx (0,x_3(t))$ and the points (0,0) and (0,1). Then we define the width of the layer as $W_{\rm s}={\rm max} \{\vert x_3\vert\}$; $W_{\rm l}={\rm min}\{\vert x_3\vert\}$ depending whether we start from the neighbourhood of the short-axis or long-axis orbit, respectively. The results are presented in Fig. 16 where we include the computed LCN, $\sigma$, (actually, $\sigma_{\rm ls}$ using the MEGNO) for the same orbits and for a motion time, $T(h)=10^4\ T_x(h)$, with T_x(h) given by (33).

The criterion to determine $h_{\rm c}(q)$, for which the layers are connected, is $W_{\rm s}(h_{\rm c})\ge W_{\rm l}(h_{\rm c})$. But W depends also on time: $W_{\rm s}$ (resp. $W_{\rm l}$) would increase (resp. decrease) with T. The results shown in Fig. 16 are for a given motion time, $T\approx 3\ 10^4\ T_x$, which appears to be large enough for our purpose. A dependence with respect to the initial conditions has to be expected too. Take, for instance, the plot for q=0.9. The layers are clearly connected at very high energies, beyond $h\approx 6$. The layer around the 1:1 resonance increases its width discontinuously except in some energy interval, , where it decreases. This may be due, perhaps, to the appearance of the 4:3 elliptic point on the p_x axis. The abrupt change in the width observed at $h\approx 4$, 4.6, 5.5 could be attributed to connections with the layers of the 4:3, 3:2 and 5:3 resonances, respectively. The LCN for the 1:1 stochastic layer decreases monotonically with h without major deviations. This reduction of the LCN with h is due to the exponential dependence of the period with h (see below).

The x-axis orbits look stable for small h, as the LCN indicate. For it turns unstable but the stochastic layer is rather thin, $W_{\rm l}$ is very close to 1 for The values of the LCN suggest that the layers could be connected at $h\approx 5.5$ while $W_{\rm s}$ reaches $W_{\rm l}$ at $h\approx 6$. Along the range , $W_{\rm s}$ and $W_{\rm l}$ are close one another. Therefore the critical value $h_{\rm c}(0.9)$ lies somewhere between (5.6, 5.9).

A similar picture is observed for q=0.8 and 0.71, but the overlap of the layers takes place for smaller h as q decreases. A simple inspection of the figures suggests that $h_{\rm c}(0.8)\approx 5$ and $h_{\rm c}(0.71)\approx 4$. In the last case (q=0.71) we observe a zone beyond $h_{\rm c}$, , where the layers appear unconnected. It may occur that some KAM curve actually exists, acting as a barrier to the diffusion or, perhaps, that a larger motion time would be necessary. In this direction, a few surfaces of section performed for these values of h and q, but for much larger motion times, do not provide any definitive answer. Anyway, in this interval, $W_{\rm s}$ and $W_{\rm l}$ are very similar. The main difference between the three figures is that, while for q=0.9 $W_{\rm s}$ reaches values close to 0.9, for q=0.8 and 0.71 $W_{\rm s}$ does not exceed 0.6.

An interesting point is that $\log\sigma$ for the gross stochastic layer decreases linearly with h. We fit by least squares the function

$$% \log\sigma(h,q)=c(q)h+d(q),$$ (36)

obtaining $c(0.9)\approx-0.457$, $c(0.8)\approx-0.438$, $c(0.71)\approx-0.421$; $d(0.9)\approx-0.09$, $d(0.8)\approx-0.08$, $d(0.71)\approx-0.06$. These estimations reveal that the LCN for the layer does not depend strongly on q, being the latter number sensitive to just one parameter, h. Nevertheless, in (36) we do not take into account the dependence of T_x on h. That is, T_x(h) was considered above only to make the computation of the LCN over similar time-scales for different energy surfaces. Therefore to obtain a meaningful result we should calculate $\sigma (h,q)T_x(h)$. Assuming that the values -0.44, -0.08 are representative for c(q) and d(q), respectively, then using (33) for T_x(h) and scaling in (36) $h\to h/p_0^2-\ln r_{\rm c}$ with p₀²=2 and $r_{\rm c}=0.1$, it is straightforward to get

$$\begin{eqnarray*}\sigma(h,q)T_x(h)\approx 0.3, \end{eqnarray*}$$

independent of h. Thus, for any energy level, $T_{\rm L}(h)\approx 3\ T_x(h)$, which is consistent with the values obtained in Sect. 4.2. This result is shown in Fig. 17 where we observe that $\sigma T_x$ is nearly constant, the smaller is q the larger is $\sigma T_x$. In any case, the LCN is bounded by .

Figure 17: Dependence of $\sigma (h,q)T_x(h)$ on h for the 1:1 stochastic layer and three different values of q. For the computation we use the values of the LCN presented in Fig. 16 and Tx given by (33). See text

The results presented in this section suggest that the motion in the logarithmic potential is mainly regular. As the energy increases, irregular motion appears confined to comparatively narrow stochastic layers. Connections among them occur at moderate-to-large energies and through thin filaments. An interesting fact is that the Lyapunov time, for the main stochastic zone, is rather short and almost independent of h and q. Note that, a similar calculation to that performed in Fig. 9 but for would lead, for "boxes'', to a completely chaotic scenario. But this is due to the choice of the one dimensional initial condition space on the xp_x-plane (the p_x axis). As mentioned in Sect. 4.2, initial conditions taken along a maximum circle on the 2D sphere (see Fig. 6) would provide a more realistic picture.

Another point that should be mentioned is that for large energies, say , the dynamics becomes almost independent of h, being q the relevant parameter. In this direction, Fig. 15 is representative of the general structure of the phase space associated to the logarithmic potential at different q levels. In fact, as h increases, the potential model approaches to the singular logarithmic potential, where the relevant parameter is only q (see Appendix).