It is true that the double pendulum exhibits integrable behavior, when the initial angles are very small, however, in general, it is very difficult to characterize the chaotic behavior of the double pendulum in terms of the initial angles. There are other representations which provide a clearer picture of its chaotic behavior.
The introductory section of the following articleby O. A. Richter, and reference therein describe the main characteristics of the double pendulum integrability (see here for the official journal version). I'll summarize here the main facts:
(Remarks: The numerical values correspond to a standard double pendulum of unit masses and unit rod lengths and unit acceleration due to gravity)
The total energy of the double pendulum is a constant of motion. The double pendulum possesses 4 equilibrium points corresponding to total energies $E= 0, 2, 4, 6$. The total energy determines the topology of the energy hypersurfaces, for $E<2$, the energy hypersurfaces are three spheres, while for $E>6$ , the energy hypersurfaces are three tori.
To further understand this point, in the case of very low energies, the system can be approximated (linearized) to an isotropic harmonic oscillator. The energy hypersurfaces are then of the form $x_1^2+x_2^2+p_1^2+p_2^2=E = const.$, while for the case of very large energies, the kinetic energy dominates and we can neglect gravity. In this case, there are two types of solutions of the equation of motion, one in which the outer rod rotates and the inner rod oscillates, and the second in which both rods rotate. The transition between the two types of solutions is determined by the value of the total angular $L$ momentum which becomes a constant of motion (due to the lack of gravity). For the standard double pendulum the transition occurs at $ L^2 = 2E$.
Both limits (small and large energies) correspond to integrable systems. This is well known, but here is a short explanation.
To see that the isotropic harmonic oscillator is integrable, one needs to solve the equations of motion in polar coordinates.
The polar angle just rotates with a constant angular velocity, and the radial coordinates oscillates in such a way that the trajectory has the shape of a two trip course through the donut hole drawn on the surface of a two-torus. This is the Liouville-Arnold torus (whose existence indicates the system's integrability) with respect to which the three sphere energy hypersurface is foliated.
In the high energy limit, a similar Liouville-Arnold torus exists when the inner rod oscillates and a torus generated by the two polar angles when the two angles rotate.
(Here the exact solution is more difficult, see for example, the following article by Enolskii, Pronine, and Richter.
Now, since both limits of vanishing and very high energies correspond to integrable systems, the total energy also controls the system characteristics, but the dependence here is much more complicated.
The transition from integrability to chaos and back as the total energy decreases from infinity to zero is described in figure 2 of Richter's article. There are a lot of details, but here are the main features: The figure correspond to the projection of the trajectories
onto the plane spanned by the outer rod angle and the total angular momenta. For very large angular momentum, the projections of the trajectories are horizontal lines with constant angular momenta (which is a constant of motion). As the energy is reduced, two disjoint chaotic regions are formed, the integrable trajectories correspond to rational and irrational tori, together with stable resonances. At about E = 10.352 which corresponds to the golden winding ratio, all irrational tori vanish and a transition to global chaos occurs.
The stable resonances also vanish eventually, at the low energies the resonances corresponding to the second integrable region start to appear.
You are making to the following mistakes:
Your initial displacement $d_0$ is 5 % of the variation of the respective variable, while it should be orders of magnitude lower. The Lyapunov exponent is defined for the limit of infinitesimal displacements, i.e., $d_0→0$.
Your observation time is far too short. You are looking at four oscillations of your dynamics, when you could have thousands. The Lyapunov exponent is defined as the average over the whole attractor/trajectory.
You do not consider multiple initial displacements or rescale the displacement vector. If you do not do this, the displacement will eventually grow to the size of the attractor and become meaningless.
In your specific example, the following happens:
You seem to be using a pendulum without small-angle approximation and thus the frequency of your pendulum depends on the intial condition.
Thus, your displaced system (red) has a slightly lower frequency than the original one and thus both systems dephase over time (at least that’s all they do in your short observation time). You see this as an increase of the distance which leads to your positive Lyapunov exponent.
If you had used an infinitesimal excitation, the different frequencies would had little impact and you would have to be very precise to see them. For methods calculating the Lyapunov exponent in tangent space (e.g., Benettin et al.’s), this would have been no problem at all.
If you had averaged over multiple initial displacements and, you would also have had displacements in which the phases would first align rather than diverge and I would intuitively say that on average you would have obtained a zero Lyapunov exponent (I cannot prove this though and it’s a rather meaningless excercise given the other mistakes).
I thought the single pendulum would produce negative exponents and the double pendulum would produce an exponent in the magnitude of 10 to 20.
The single pendulum should have a zero (largest) Lyapunov exponent like all limit cycle dynamics. You can understand this as follows: A displacement along the trajectory will neither shrink nor grow in time (on average). The two time series will be exactly the same except for an offset corresponding to the displacement.
I thought […] the double pendulum would produce an exponent in the magnitude of 10 to 20.
The absolute value of the Lyapunov exponent is rather meaningless, as it depends on how you scale time. Without knowing the detailed parameters of your double pendulum, it’s impossible to say what you should expect. The only thing that is meaningful is the sign of the largest Lyapunov exponent.
Best Answer
I found the (shamefully simple) answer myself:
The simple double pendulum is a conservative system, hence due to Liouville’s theorem the phase-space volume given by a given ensemble of trajectories is constant over time. However, if the system exhibited transients and thus attractors, all trajectories starting within the basin of attraction of a given attractor would eventually converge towards that attractor. As the phase-space volume of the basin of attraction is larger than that of the attractor, Liouville’s theorem would be broken. Thus transients and attractors can only occur in dissipative systems.