In the study of Hamiltonian systems finding conserved quantities may give precious insights about the dynamics of the system.   One of this quantity is the so-called adiabatic invariant, which is essentially the area subtended by one periodic trajectory in the canonical phase space.  In our case, the adiabatic invariant would be:

(where G=sqrt((1-e^2)*(\mu*a)) and g=\omega) over one trajectory.    In our case, we computed the area subtended by the separatrix for several values of H=sqrt((1-e^2)*(\mu*a))*cos(i).  The next figure show the separatrix for H/(sqrt(/mu*a)=0.0030.

We computed the area covered by the separatrix for different values of H.  Results are in the following figures: in the first we report the dependence of the fraction of phase space covered by the libration island with respect to the whole available phase space (both libration islands at + and -90o are considered), with respect to H/(sqrt(/mu*a); in the second the value of J/(sqrt(/mu*a) with respect to H/(sqrt(/mu*a) (one libration island only), and in the third J vs. H (sqrt(/mu*a)=2.06982*10¹³ m²/s, for S/2000S5 (a=0.075 AU, again only one libration island)).
 

In the next image we show a plot of the adiabatic invarant versus H for different values of a.

In the next figures we report the plots of J vs. H for a numerical simulation involving S5 (no drag involved).  The plots on the left are J vs. H, J vs. a, and a vs. H, respectively.
 

In the next two plots we report the time evolution of orbital elements for  a simulation involving S5, with Stokes drag, and J vs. H for the same simulation (J was computed until the time the particle left the Kozai resonance.  The colored line report the values of J vs. H for the minimum, average and maximum value of a, respectively.