About Tidal Evolution of Quasi-Periodic Orbits of Satellites

Abstract

Tidal interactions between Planet and its satellites are known to be the main phenomena, which are determining the orbital evolution of the satellites. The modern ansatz in the theory of tidal dissipation in Saturn was developed previously by the international team of scientists from various countries in the field of celestial mechanics. Our applying to the theory of tidal dissipation concerns the investigating of the system of ODE-equations (ordinary differential equations) that govern the orbital evolution of the satellites; such an extremely non-linear system of 2 ordinary differential equations describes the mutual internal dynamics for the eccentricity of the orbit along with involving the semi-major axis of the proper satellite into such a monstrous equations. In our derivation, we have presented the elegant analytical solutions to the system above; so, the motivation of our ansatz is to transform the previously presented system of equations to the convenient form, in which the minimum of numerical calculations are required to obtain the final solutions. Preferably, it should be the analytical solutions; we have presented the solution as a set of quasi-periodic cycles via re-inversing of the proper ultra-elliptical integral. It means a quasi-periodic character of the evolution of the eccentricity, of the semi-major axis for the satellite orbit as well as of the quasi-periodic character of the tidal dissipation in the Planet.

1 Introduction, the System of Equations

Recently, an epochal analyses were executed over a lot of data of circa 130-years observations regarding the Saturn positions in a space, accomplished with the very accurate observations of positions for all the satellites of Saturn. It was a hard work to find a proper information (in various sources of data) about the appropriate observations, to check their validity, then to combine it for future analysis in the united data-base for computations. Such an analyses were made by the international group of scientists in the comprehensive articles (Lainey et al. 2012, 2015). Authors used a numerical methods to obtain a theoretical solutions, then they compared it with the data of all the reasonable observations. For example, in Appendix of Lainey et al. (2012, 2015) a systems of equations (A1), (A2) have been stated for mutual evolution of the eccentricity e along with the semi-major axis a of the moons of Saturn. Here and below we note that the tidal effects are introduced by means of the Love number k₂, which is describing the response of the potential of the distorted body in regard to the experiencing tides, as well as by the quality factor Q, which is inversely proportional to the amount of energy dissipated essentially as heat by tidal friction (Lainey et al. 2009); so, tidal effects are introduced in the combination k₂/Q for Saturn and satellite.

In particular, we recall that we have (as a first approximation) for the tides raised in the primary (a case of tidal interaction between Saturn and Titan in Kaula (1964) was considered):

$$\begin{aligned} \frac{da}{dt} & {\kern 1pt} = \frac{{3k_{2} mnR^{5} }}{{QMa^{4} }}\left( {1 + \frac{51}{4}e^{2} } \right), \\ \frac{de}{dt} & {\kern 1pt} = \frac{{57k_{2} mn}}{8QM}\left( {\frac{R}{a}} \right)^{5} e, \\ \end{aligned}$$

(A1)

here m is the mass of satellite, M is the mass of Saturn, n is the osculating mean motion, R is equatorial radius.

But we also recall that we have (as a first approximation) for the tides raised in the 1:1 spin–orbit satellite (Peale and Cassen 1978):

$$\begin{aligned} \frac{da}{dt}{\kern 1pt} & = - \frac{{21k_{2}^{s} MnR_{s}^{5} }}{{Q^{s} ma^{4} }}e^{2} , \\ \frac{de}{dt} & = - \frac{{21k_{2}^{s} Mn}}{{2Q^{s} m}}\left( {\frac{{R_{s} }}{a}} \right)^{5} e \\ \end{aligned}$$

(A2)

here sign “s” denotes the case of satellite.

Besides, we should note that, according to the Kepler’s law of orbital motion, the square of mean motion:

$$n^{2} = \frac{G(M + m)}{{a^{3} }},$$

(*)

where G—gravitational constant of the Newton’s law of universal gravitation.

2 Reduction of the System of Equations (A2)

We could present system of equations (A2) as below

$$\begin{aligned} \frac{da}{dt}{\kern 1pt} & = - \frac{{21k_{2}^{s} M\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{{\kern 1pt} \frac{3}{2}}} }}} \right) \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}}\left( {\frac{1}{a}} \right)^{4} \cdot e^{2} , \\ \frac{de}{dt} & = - \frac{{21k_{2}^{s} M\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{\frac{3}{2}}} }}} \right) \cdot \left( {R_{s} } \right)^{5} }}{{2Q^{s} m}}\left( {\frac{1}{a}} \right)^{5} {\kern 1pt} \cdot e \\ \end{aligned}$$

(2.1)

Let us denote just for simplicity

$$B = \frac{{21k_{2}^{s} M\left( { \pm \sqrt {G(M + m)} } \right) \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}}.$$

Mathematical procedure for reduction of the system of Eqs. (2.1) has been moved to an “Appendix”, with only the resulting formulae left in the main text ({a ₀, e ₀} = {a(0), e(0)} = const):

$$a{\kern 1pt} = a_{0} \cdot \exp (e^{2} - e_{0}^{2} ),$$

(2.2)

Thus, if we consider the case of eccentricity e → 0, we could obtain in “Appendix” from Eqs. (2.3)–(2.5) (here below Δt should be considered as long time-period scale):

$$e \cong e_{0} \cdot \exp \left( { - \frac{B}{2} \cdot \frac{{\exp \left( {\frac{13}{2}e_{0}^{2} } \right)}}{{(a_{{{\kern 1pt} 0}} )^{{\frac{13}{2}}} }} \cdot \Delta t} \right),$$

(2.6)

just compare it with the appropriate plot at Fig. 3 in Lainey et al. (2012); as well as we could obtain the appropriate expression for the semi-major axis from (2.2), using (2.6):

$$a{\kern 1pt} = a_{0} \cdot \exp \left( {e_{0} \cdot \exp \left( { - \frac{B}{2} \cdot \frac{{\exp \left( {\frac{13}{2}e_{0}^{2} } \right)}}{{(a_{0} )^{{\frac{13}{2}}} }} \cdot \Delta t} \right)} \right)$$

(2.7)

where the scale-factor a ₀ should be given by the initial conditions. We schematically imagine the plot of solution (2.7) at Fig. 1 as presented below.

Fig. 1

Schematically imagined the plot of solution (2.7) for the function a(t)

3 Reduction of the System of Equations (A1)

Let us consider the more complicated case of system of equations (A1) as below

$$\begin{aligned} \frac{da}{dt}{\kern 1pt} & = \frac{{3k_{2} m\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{{\kern 1pt} \frac{3}{2}}} }}} \right)R^{5} }}{{QMa^{4} }}\left( {1 + \frac{51}{4}e^{2} } \right), \\ \frac{de}{dt} & = \frac{{57k_{2} m\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{\frac{3}{2}}} }}} \right)R^{5} }}{8QM}\left( {\frac{1}{a}} \right)^{5} \cdot e, \\ \end{aligned}$$

(3.1)

Let us also denote just for simplicity

$$A = \frac{{k_{2} m\left( { \pm \sqrt {G(M + m)} } \right) \cdot R^{5} }}{QM}.$$

Mathematical procedure for reduction of the system of Eq. (3.1) has also been moved to an “Appendix”, with only the resulting formulae left in the main text ({a ₁ , e ₁} = {a(0), e(0)} = const):

$$a{\kern 1pt} = a_{1} \cdot {\kern 1pt} \left( {\frac{e}{{e_{1} }}} \right)^{{\frac{8}{19}}} \cdot \exp \left( {\frac{51}{19}(e^{{{\kern 1pt} 2}} - e_{1}^{2} )} \right),$$

(3.2)

where the term: exp ((51/19)·(e ² − e ²₁ )) ≅ 1. So, using (3.2), we could obtain from the 2nd of Eqs. (3.1) (here below Δt should be considered as long time–period scale):

$$e = e_{1} \cdot \left( {\frac{39}{2} \cdot A \cdot \left( {a_{1} } \right)^{{ - {\kern 1pt} \frac{13}{2}}} \cdot \Delta t} \right)^{{\frac{19}{52}}}$$

(3.3)

where the scale-factor a ₁ should be given by the initial conditions according to the assumption e → 0.

We schematically imagine an approximation of the solution (3.3) dynamics at Fig. 2 (where we assume the extent (19/52) ≅ 1/2 just for simplicity of presentation below):

Fig. 2

Schematically imagined the plot of solution (3.3) for the function e(t)

Analyzing the contributions (2.6) and (3.3) into the effect of tidal dissipation (in regard to the evolution of the orbit of satellite motion around the planet), we can see that such the contributions apparently differ from each other.

Indeed, the contribution of tidal dissipation in planet (2.6)–(2.7) tends to decrease the eccentricity as well as the semi-major axis of the satellite orbit [depending on the sign of mean motion (*)], but a proper contribution of tidal dissipation in satellite (3.2)–(3.3) tends to increase the eccentricity as well as the semi-major axis of the satellite orbit [and vice versa, depending on the sign of mean motion (*)].

Thus, we should evaluate the combined contribution of the effects of tidal dissipation both in the planet and satellite.

4 Reduction of the Combined System of Equations (A1) + (A2)

Let us consider the case of the combined system of equations (A1) + (A2) in the sense of combined contributions to the tidal dissipation (of the planet + satellite) as below

$$\begin{aligned} \frac{da}{dt}{\kern 1pt} & = \frac{{3k_{2} m\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{\frac{3}{2}}} }}} \right)R^{5} }}{{QMa^{{{\kern 1pt} 4}} }}\left( {1 + \frac{51}{4}e^{2} } \right) - \frac{{21k_{2}^{s} M\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{\frac{3}{2}}} }}} \right) \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}}\left( {\frac{1}{a}} \right)^{4} {\kern 1pt} \cdot e^{2} , \\ \frac{de}{dt} & = \frac{{57k_{2} m\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{\frac{3}{2}}} }}} \right)R^{5} }}{8QM}\left( {\frac{1}{a}} \right)^{\,5} \cdot e - \frac{{21k_{2}^{s} M\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{\frac{3}{2}}} }}} \right) \cdot \left( {R_{s} } \right)^{5} }}{{2Q^{s} m}}\left( {\frac{1}{a}} \right)^{5} {\kern 1pt} \cdot e. \\ \end{aligned}$$

The last system could be transformed to the form below

$$\begin{aligned} \frac{da}{dt} & = \left( {\frac{C}{{a^{{\frac{11}{2}}} }}} \right) \cdot \left( {D + E \cdot e^{2} } \right), \\ \frac{de}{dt} & = \left( {\frac{C}{{2a^{{\frac{13}{2}}} }}} \right) \cdot F \cdot e, \\ \end{aligned}$$

(4.1)

where we have denoted (just for simplicity) the appropriate constants:

$$\begin{aligned} C & = \pm 3\sqrt {G(M + m)} {\kern 1pt} ,\quad D = \frac{{k_{2} m \cdot R^{5} }}{QM},\quad E = \left( {\frac{{k_{2} m \cdot R^{5} }}{QM} \cdot \frac{51}{4} - \frac{{7k_{2}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}}} \right), \\ F & = \left( {\frac{{19k_{2} mR^{5} }}{4QM} - \frac{{7k_{2}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}}} \right). \\ \end{aligned}$$

(4.2)

Also, the mathematical procedure for reduction of the system (3.1) has been moved to an “Appendix”, with only the resulting formulae left in the main text ({a ₂ , e ₂} = {a(0), e(0)} = const):

$$a = a_{2} \cdot {\kern 1pt} \left( {\frac{e}{{e_{2} }}} \right)^{{2\,\left( {\frac{D}{F}} \right)}} \cdot \exp \left( {\left( {\frac{E}{F}} \right) \cdot (e^{2} - e_{2}^{2} )} \right),$$

(4.3)

where the term: exp ((E/F)·(e ² − e ²₂ )) ≅ 1. So, using (4.3), we could obtain from the 2nd of Eq. (4.1) (here below Δt should be considered as long time-period scale):

$$e = e_{2} \cdot \left( {\frac{13C \cdot D}{2} \cdot \left( {a_{2} } \right)^{{ - \frac{13}{2}}} \cdot \Delta t} \right)^{{\frac{F}{13D}}}$$

(4.4)

where the scale-factor a ₂ should be given by the initial conditions according to the assumption e → 0.

If we assume the extent (F/(13D)) ≅ 1/2 in (4.4) (just for simplicity of presentation), we should conclude that the approximate dynamics of the solution (4.4) coincide to the dynamics of previously discussed solution (3.3) which was imagined at Fig. 2.

The last but not least, we should especially note that the resulting combined contributions of tidal dissipation of both the Planet and satellite depend on the ratio (F/D) (4.2) (see the appropriate expressions for solution (4.3)–(4.4)), which is obviously not depending on the sign of mean motion (*):

$$F = \left( {\frac{{19k_{2} mR^{5} }}{4QM} - \frac{{7k_{2}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}}} \right),\quad D = \frac{{k_{2} m \cdot R^{5} }}{QM}$$

5 Discussions

Tidal interactions between Planet and its satellites are known to be the main phenomena, which are determining the orbital evolution of the satellites. There are a lot of theories of tidal dissipation, but most of them could be associated with two main types: (1) tidal friction for bodies with fluid layers, (2) solid tidal dissipation. Indeed, the rheological law for the actual rheological parameters, obeyed by the material of the bodies, and their role in dissipation differ from one aforementioned types to another.

Definitely, a short review on the applicability over formulations for solid tidal dissipation would be helpful. In this respect we confine ourselves to mention the paper (Efroimsky and Makarov 2013) in which most popular cases of tidal friction are remarked (as the constant geometric lag model or the constant time lag model).

This paper presents a mathematical technique that helps with the analytical integration over time of eccentricity + semi-major axis in a binary system experiencing tidal friction. We should note that integration of these equations is quite fundamental to many studies and so even small improvements in the method can be a benefit.

The described formulation incorporates formulae for the tidal friction that is surely not appropriate for the aforementioned bodies with fluid layers. In this respect we confine ourselves to mention the paper (Tyler 2014) in which all the difficulties concerning the most complicated cases of tidal friction for bodies with fluid layers are remarked. A leading result in the comprehensive study above is that the tidal response expected cannot simply be inferred from the orbit, or even the expected Q (quality factor). Also, the differences in the nature of the fluid versus solid tidal dissipation have been pointed out in the aforementioned article (the dependence on ocean thickness is also at least as important as Q). Referring to the comprehensive article (Tyler 2014), we should generalize our future researches for tidal dissipation effect, where we should consider or describe what systems their equations apply to and whether they can be extended to the case of bodies with fluids.

6 Conclusion

Our applying to the theory of tidal dissipation concerns the investigating of the system of ODE-equations that govern the orbital evolution of the satellites; such an extremely non-linear system of 2 ordinary differential equations describes the mutual internal dynamics for the eccentricity of the orbit along with involving the semi-major axis of the proper satellite to such a monstrous equations.

Referring to the comprehensive articles (Efroimsky and Lainey 2007; Efroimsky 2015) we should generalize our future researches for tidal dissipation effect depending on the tidal-flexure frequency χ.

Indeed, according to the modern ansatz (Efroimsky and Lainey 2007) the quality factor Q of the Planet could be assumed depending on the tidal-flexure frequency χ as below:

$$Q\, \propto \chi^{\alpha } ,\quad \alpha = 0.16 \div 0.4$$

where frequency χ is apparently depending on the mean motion (*): χ = 2 |ω_p − n|, according to results reported in Efroimsky and Lainey (2007) (ω_p being planet’s spin rate).

Besides, the ratio ($k_{2}^{s} /Q^{s}$) of the satellite could be assumed depending on the tidal-flexure frequency χ as below (Efroimsky 2015):

$$\frac{{k_{2}^{s} }}{{Q^{s} }} \propto (\eta \cdot \chi )^{\beta } ,\quad \beta = \pm 1$$

where η is the effective viscosity of the satellite; but frequency χ is apparently also depending on the mean motion (*) as above.

We should recall that initial system of Eq. (A2) was presented for the tides raised in the 1:1 resonance for the spin–orbit of satellite. So, we should up-date it for the case where frequency χ is supposed to be depending on the mean motion (*): χ = 2 |ω_p − n|, according to the ansatz Efroimsky and Lainey (2007) (ω _p being the planet’s spin rate).

It means that we should correct properly the set of coefficients {C, D, E, F} in formulae (4.1)–(4.2) for such a case; especially, the set {D, F} should be corrected as the coefficients, which are determining the structure of the solution (4.3)–(4.4) across the ratio (D/F) or (F/D) as the key dynamical parameter of the system. Meanwhile, such a correction could be accomplished with the data of astrometric observations: indeed, we could adjust analytical solutions with respect to the actual data of observations for the orbits of satellites.

All in all, the physically reasonable hypothesis should be assumed as below (we assume the mixed scenario): the tidal dissipation of the satellite is assumed to be equal to the constant value, but the tidal dissipation of the Planet could be assumed depending on the tidal-flexure frequency χ as suggested above:

$$\begin{aligned} Q \propto\upchi^{\upalpha} ,\quad (\upalpha = 0.2 \div 0.4)\quad n^{2} = \frac{G(M + m)}{{a^{3} }}, \hfill \\ \quad \Rightarrow\upchi = 2\left| {\upomega_{p} - \sqrt {\frac{G(M + m)}{{a^{3} }}} } \right| \hfill \\ \end{aligned}$$

(6.1)

In such a case, we could obtain from the Eqs. (4.1), (6.1) (see the proper derivation (6.2)–(6.9) in “Appendix”), with only the resulting formulae left in the main text below:

$$\int {\frac{{\left( {\upomega_{p} - \sqrt {\frac{G(M + m)}{{(a(0))^{3} \cdot \left( {\frac{e}{e(0)}} \right)^{{\frac{3}{H}}} }}} } \right)^{\alpha } \cdot \left( {\frac{e}{e(0)}} \right)^{{\,{\kern 1pt} \left( {\frac{13}{2H} - 1} \right)}} }}{{\left( {\frac{{19k_{2} mR^{5} }}{{4M \cdot 2^{\upalpha} }} - \frac{{7k_{2}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}} \cdot \left( {\upomega_{p} - \sqrt {\frac{G(M + m)}{{(a(0))^{3} \cdot \left( {\frac{e}{e(0)}} \right)^{{\frac{3}{H}}} }}} } \right)^{\upalpha} } \right)}}d\left( {\frac{e}{e(0)}} \right)} = \left( {\frac{C}{{2(a(0))^{{\frac{13}{2}}} }}} \right) \cdot \int {dt}$$

(6.6)

which could be simplified by Tailor series (as first approximation) in regard to the term below

$$\left( {\upomega_{p} - \sqrt {\frac{G(M + m)}{{(a(0))^{3} \cdot \left( {\frac{e}{e(0)}} \right)^{{{\kern 1pt} \frac{3}{H}}} }}} } \right)^{\upalpha} \cong \left( {\upomega_{p} } \right)^{\upalpha} \cdot \left( {1 - \left( {\frac{\alpha }{{\omega_{p} }}} \right) \cdot \sqrt {\frac{G(M + m)}{{(a(0))^{3} \cdot \left( {\frac{e}{e(0)}} \right)^{{\frac{3}{H}}} }}} } \right).$$

So, we obtain ({a(0), e(0)} = const):

$$\begin{aligned} a \cong a(0) \cdot u^{{ - \frac{2}{3}}} ,\quad e = e(0) \cdot u^{{ - \frac{2H}{3}}} , \hfill \\ H = \left( {\frac{19}{8} - \frac{{7k_{2}^{s} M^{2} \cdot \left( {R_{s} } \right)^{5} }}{{2k_{2} \cdot Q^{s} m^{2} R^{5} }} \cdot \left( {2\omega_{p} } \right)^{\alpha } } \right), \hfill \\ \end{aligned}$$

$$\begin{aligned} \int {\frac{{u^{{{\kern 1pt} \left({- \frac{16}{3}} \right)}}}}{{\left({\varLambda \cdot \left({1 \pm \left({\frac{\alpha C}{{3\omega_{p} \cdot (a(0))^{{\frac{3}{2}}}}}} \right) \cdot u} \right) - \varXi} \right)}}du = - \left({\frac{3C}{{2H \cdot (a(0))^{{\frac{13}{2}}}}}} \right) \cdot \int {dt},} \hfill \\ \varLambda = \frac{{19k_{2} mR^{5}}}{{8\left({2\omega_{p}} \right)^{\alpha} \cdot M}},\quad C = \pm 3\sqrt {G(M + m)},\quad \varXi = \frac{{7k_{2}^{s} M \cdot \left({R_{s}} \right)^{5}}}{{2Q^{s} m}} \hfill \\ \end{aligned}$$

(6.9)

Solution (6.9) could obviously be reduced to the solution (4.4) if we choose α = 0; in this case we should note that

$$\left( {2\upomega_{{{\kern 1pt} p}} } \right)^{\upalpha} \to Q = const.$$

Let us note that the left part of Eq. (6.6) is the proper ultra-elliptical integral of fractional order in regard to the function u (see Lawden 1989), which depends on the eccentricity e in formulae (6.9). But the elliptical integral is known to be a generalization of a class of inverse periodic functions. Thus, by the obtaining of re-inverse dependence for the expression (6.6), we could present the solution as a set of quasi-periodic cycles: it means a quasi-periodic character of the evolution of the eccentricity, of the semi-major axis for the satellite orbit as well as of the quasi-periodic character of the tidal dissipation in the Planet.

By the way, if we take into consideration the dependence of the ratio ($k_{2}^{s} /Q^{s}$ k ₂ˢ/Qˢ) for the satellite on the tidal-flexure frequency χ as below (Efroimsky 2015):

$$\frac{{k_{2}^{s} }}{{Q^{s} }} \propto (\upeta \cdot\upchi)^{\upbeta} ,\quad\upbeta = \pm 1$$

recall that χ = |n| for the case of satellite (Lainey et al. 2012, 2015), we should also obtain the ultra-elliptical integral in regard to the function a as the analogue of Eq. (6.2) for such a case. Nevertheless, we restrict ourselves to the chosen case of constant tidal dissipation inside the satellite for the current research (the satellite is assumed synchronised).

Finally, we should especially note that the mean motion (*) is evaluated here according to the Kepler’s law of orbital motion. But in the case of restricted 3-bodies problem (RTBP) Ershkov (2015, 2017) it should differ from the case of classical solution of two-bodies problem.

Appendices

Appendix (Calculations in 2–6 Sections)

1.1 Reduction of the System of Equations (A2)

We could present system of equations (A2) as below

$$\begin{aligned} \frac{da}{dt} & {\kern 1pt} = - \frac{{21k_{2}^{s} M\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{\frac{3}{2}}} }}} \right) \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}}\left( {\frac{1}{a}} \right)^{4} \cdot e^{2} , \\ \frac{de}{dt} & = - \frac{{21k_{2}^{s} M\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{{\kern 1pt} \frac{3}{2}}} }}} \right) \cdot \left( {R_{s} } \right)^{5} }}{{2Q^{s} m}}\left( {\frac{1}{a}} \right)^{5} {\kern 1pt} \cdot e \\ \end{aligned}$$

(2.1)

Let us denote just for simplicity

$$B = \frac{{21k_{2}^{s} M\left( { \pm \sqrt {G(M + m)} } \right) \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}},$$

if we divide each part of 1st equation on the proper part of the 2nd equation of (2.1), it should yield ({a ₀, e ₀} = {a(0), e(0)} = const):

$$\begin{aligned} \frac{da}{de}{\kern 1pt} = 2\frac{e}{(1/a)}, \Rightarrow \frac{da}{a} = d(e^{2} ), \Rightarrow \ln \,a{\kern 1pt} = e^{2} + const, \Rightarrow a = a_{0} \cdot \exp \,(e^{2} - e_{0}^{2} ), \hfill \\ \end{aligned}$$

(2.2)

so, using (2.2), we could obtain from the 2nd of Eqs. (2.1)

$$\begin{aligned} \frac{de}{dt} = - \left( {\frac{B}{2}} \right) \cdot e \cdot a^{{ - \frac{13}{2}}} \Rightarrow \frac{de}{dt} = - \left( {\frac{B}{2}} \right) \cdot e{\kern 1pt} \cdot (a_{0} \cdot \exp (e^{2} - e_{0}^{2} ))^{{ - \frac{13}{2}}} \hfill \\ \Rightarrow \int {\left( {\frac{{\exp \left( {\frac{13}{2}e^{2} } \right)}}{e}} \right){\kern 1pt} de} = - \left( {\frac{B}{2}} \right) \cdot \frac{{\exp \left( {\frac{13}{2}e_{0}^{2} } \right)}}{{(a_{0} )^{{\frac{13}{2}}} }} \cdot \int {dt} \hfill \\ \end{aligned}$$

(2.3)

The left side of expression (2.3) could be transformed by the appropriate change of variables as below

$$\begin{aligned} e{\kern 1pt} = \sqrt u , \Rightarrow \int {\left( {\frac{{\exp \left( {\frac{13}{2}u} \right)}}{{2\left( {\sqrt u } \right)^{2} }}} \right)} du{\kern 1pt} = - \left( {\frac{B}{2}} \right) \cdot \frac{{\exp \left( {\frac{13}{2}e_{0}^{2} } \right)}}{{(a_{{{\kern 1pt} 0}} )^{{\frac{13}{2}}} }} \cdot \int {dt} \hfill \\ \Rightarrow \int {\left( {\frac{{\exp \left( {\frac{13}{2}u} \right)}}{u}} \right)} du = - B \cdot \frac{{\exp \left( {\frac{13}{2}e_{0}^{2} } \right)}}{{(a_{0} )^{{\frac{13}{2}}} }}\int {dt} \hfill \\ \end{aligned}$$

(2.4)

where (Kamke 1971)

$$\int {\left( {\frac{{\exp \left( {\frac{13}{2}u} \right)}}{u}} \right)} du = \ln \,\left| u \right| + \frac{{\left( {\frac{13}{2}} \right) \cdot u}}{1 \cdot 1!} + \frac{{\left( {\frac{13}{2}} \right)^{2} \cdot u^{2} }}{2 \cdot 2!} + \frac{{\left( {\frac{13}{2}} \right)^{3} \cdot u^{3} }}{3 \cdot 3!} + \cdots + \frac{{\left( {\frac{13}{2}} \right)^{n} \cdot u^{n} }}{n \cdot n!} + \ldots ,\quad u = e^{2}$$

(2.5)

Thus, if we consider the case of eccentricity e → 0, we could obtain from Eqs. (2.4)–(2.5):

$$2\ln \left( {\frac{e}{{e_{0} }}} \right) \cong - B \cdot \frac{{\exp \left( {\frac{13}{2}e_{0}^{2} } \right)}}{{(a_{0} )^{{\frac{13}{2}}} }} \cdot \Delta t \Rightarrow e \cong e_{0} \cdot \exp \left( { - \frac{B}{2} \cdot \frac{{\exp \left( {\frac{13}{2}e_{{{\kern 1pt} 0}}^{2} } \right)}}{{(a_{{{\kern 1pt} 0}} )^{{\frac{13}{2}}} }} \cdot \Delta t} \right),$$

(2.6)

just compare it with the appropriate plot at Fig. 3 in Lainey et al. (2012); as well as we could obtain from Eqs. (2.4)–(2.5) the appropriate expression for the semi-major axis from (2.2) (here below Δt should be considered as long time-period scale):

$$a{\kern 1pt} = a_{0} {\kern 1pt} \cdot \exp \left( {e_{0} \cdot \exp \left( { - \frac{B}{2} \cdot \frac{{\exp \left( {\frac{13}{2}e_{0}^{2} } \right)}}{{(a_{0} )^{{\frac{13}{2}}} }} \cdot \Delta t} \right)} \right)$$

(2.7)

where the scale-factor a ₀ should be given by the initial conditions.

1.2 Reduction of the System of Equations (A1)

Let us consider the more complicated case of system of equations (A1) as below

$$\begin{aligned} \frac{da}{dt}{\kern 1pt} = \frac{{3k_{2} m\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{{\kern 1pt} \frac{3}{2}}} }}} \right)R^{5} }}{{Q\,M\,a^{{{\kern 1pt} 4}} }}\left( {1 + \frac{51}{4}e^{2} } \right), \hfill \\ \frac{de}{dt}{\kern 1pt} = \frac{{57k_{2} m\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{{\kern 1pt} \frac{3}{2}}} }}} \right)R^{5} }}{8QM}\left( {\frac{1}{a}} \right)^{5} \cdot e, \hfill \\ \end{aligned}$$

(3.1)

Let us also denote just for simplicity

$$A\;\; = \;\;\frac{{k_{{{\kern 1pt} 2}}^{\kern 1pt} \,m\,\left( { \pm \sqrt {G\,(M + m)} } \right){\kern 1pt} \, \cdot \,R_{\,}^{5} }}{Q\,\,M\,}\;,$$

then, if we divide each part of 1st equation on the proper part of 2nd equation of (3.1), it should yield ({a ₁ , e ₁} = {a(0), e(0)} = const):

$$\begin{aligned} \frac{da}{de}{\kern 1pt} = \frac{8}{19}\frac{{\left( {1 + \frac{51}{4}e^{2} } \right)}}{(1/a) \cdot e}, \Rightarrow \frac{da}{a} = \frac{8}{19}\frac{{\left( {1 + \frac{51}{4}e^{2} } \right)}}{e}de, \hfill \\ \ln \,a{\kern 1pt} = \frac{8}{19}\ln \,e + \frac{51}{19}e^{2} + const, \Rightarrow a = a_{1} \cdot {\kern 1pt} \left( {\frac{e}{{e_{1} }}} \right)^{{\frac{8}{19}}} \cdot \exp \left( {\frac{51}{19}(e^{{{\kern 1pt} 2}} - e_{1}^{2} )} \right), \hfill \\ \end{aligned}$$

(3.2)

where the term: exp ((51/19)·(e ² − e ²₁ )) ≅ 1. So, using (3.2), we could obtain from the 2nd of Eqs. (3.1) (here below Δt should be considered as long time-period scale):

$$\begin{aligned} \frac{{d\left( {\frac{e}{{e_{1} }}} \right)}}{dt} = \frac{57}{8}A\left( {a_{1} \cdot \left( {\left( {\frac{e}{{e_{1} }}} \right)^{{\frac{8}{19}}} } \right)} \right)^{{ - \frac{13}{2}}} \cdot \left( {\frac{e}{{e_{1} }}} \right), \Rightarrow \int {\left( {\frac{e}{{e_{1} }}} \right)^{{\frac{33}{19}}} } d\left( {\frac{e}{{e_{1} }}} \right) = \frac{57}{8}A \cdot \left( {a_{1} } \right)^{{ - \frac{13}{2}}} \cdot \int {dt} , \Rightarrow \hfill \\ \frac{1}{{\left( {\frac{33}{19} + 1} \right)}}\left( {\frac{e}{{e_{1} }}} \right)^{{\frac{52}{19}}} = \frac{57}{8} \cdot A \cdot \left( {a_{1} } \right)^{{ - \frac{13}{2}}} \cdot \int {dt} \Rightarrow e = e_{1} \cdot \left( {\frac{39}{2} \cdot A \cdot \left( {a_{1} } \right)^{{ - \frac{13}{2}}} \cdot \Delta t} \right)^{{\frac{19}{52}}} \hfill \\ \end{aligned}$$

(3.3)

where the scale-factor a ₁ should be given by the initial conditions according to the assumption e → 0.

1.3 Reduction of the Combined System of Equations (A1) + (A2)

Let us consider the case of the combined system of equations (A1) + (A2) in the sense of combined contributions to the tidal dissipation (of Saturn + satellite) as below

$$\begin{aligned} \frac{da}{dt}{\kern 1pt} & = \frac{{3k_{2} m\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{\frac{3}{2}}} }}} \right)R^{5} }}{{QMa^{4} }}\left( {1 + \frac{51}{4}e^{2} } \right) - \frac{{21k_{2}^{s} M\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{\frac{3}{2}}} }}} \right) \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} \,m}}\left( {\frac{1}{a}} \right)^{4} {\kern 1pt} \cdot e^{2} , \\ \frac{de}{dt} & = \frac{{57k_{2} m\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{{\kern 1pt} \frac{3}{2}}} }}} \right)R^{5} }}{8QM}\left( {\frac{1}{a}} \right)^{5} \cdot e - \frac{{21k_{2}^{s} M\left( { \pm \frac{{\sqrt {G(M + m)} }}{{a^{{{\kern 1pt} \frac{3}{2}}} }}} \right) \cdot \left( {R_{s} } \right)^{5} }}{{2Q^{s} m}}\left( {\frac{1}{a}} \right)^{5} {\kern 1pt} \cdot e. \\ \end{aligned}$$

The last system could be transformed to the form below

$$\begin{aligned} \frac{da}{dt}{\kern 1pt} & = \left( {\frac{C}{{a^{{\frac{11}{2}}} }}} \right) \cdot \left( {D + E \cdot e^{2} } \right), \\ \frac{de}{dt} & = \left( {\frac{C}{{2a^{{\frac{13}{2}}} }}} \right) \cdot F \cdot e, \\ \end{aligned}$$

(4.1)

where we have denoted (just for simplicity) the appropiate constants:

$$\begin{aligned} C{\kern 1pt} & = \pm 3\sqrt {G(M + m)},\quad D = \frac{{k_{2} m \cdot R^{5} }}{QM}, \quad E{\kern 1pt} = \left( {\frac{{k_{2} m \cdot R^{5} }}{QM} \cdot \frac{51}{4} - \frac{{7k_{{{\kern 1pt} 2}}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}}} \right), \\ F & = \left( {\frac{{19k_{2} mR^{5} }}{4QM} - \frac{{7k_{2}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}}} \right). \\ \end{aligned}$$

(4.2)

Then, if we divide each part of 1st equation on the proper part of 2nd equation of (4.1), it should yield ({a ₂ , e ₂} = {a(0), e(0)} = const):

$$\begin{aligned} \frac{da}{de}{\kern 1pt} & = 2\frac{{\left( {D + E \cdot e^{2} } \right)}}{(1/a) \cdot F \cdot e}, \Rightarrow \frac{da}{a} = \left( {2\left( {\frac{D}{F}} \right) \cdot \frac{1}{e} + 2\left( {\frac{E}{F}} \right) \cdot e} \right)de, \\ \ln \,a & = 2\left( {\frac{D}{F}} \right)\,\ln \,e + \left( {\frac{E}{F}} \right) \cdot e^{2} + const, \Rightarrow a = a_{2} \cdot \left( {\frac{e}{{e_{2} }}} \right)^{{2\,\left( {\frac{D}{F}} \right)}} \cdot \exp \left( {\left( {\frac{E}{F}} \right) \cdot (e^{2} - e_{2}^{2} )} \right), \\ \end{aligned}$$

(4.3)

where the term: exp ((E/F)·(e ² − e ²₂ )) ≅ 1. So, using (4.3), we could obtain from the 2nd of Eqs. (4.1) (here below Δt should be considered as long time-period scale):

$$\begin{aligned} \frac{{d{\kern 1pt} \left( {\frac{e}{{e_{2} }}} \right)}}{dt} & = \frac{C}{{2\left( {a_{2} \cdot {\kern 1pt} \left( {\frac{e}{{e_{{{\kern 1pt} 2}} }}} \right)^{{2\,\left( {\frac{D}{F}} \right)}} } \right)^{{\frac{13}{2}}} }} \cdot F \cdot \left( {\frac{e}{{e_{2} }}} \right), \Rightarrow \int {\left( {{\kern 1pt} \left( {\frac{e}{{e_{{{\kern 1pt} 2}} }}} \right)^{{\left( {13\,\left( {\frac{D}{F}} \right) - 1} \right)}} } \right)} d{\kern 1pt} \left( {\frac{e}{{e_{2} }}} \right) = \frac{C}{2} \cdot F \cdot \left( {a_{2} } \right)^{{ - \frac{13}{2}}} \cdot \int {dt} , \Rightarrow \\ \frac{1}{{13\left( {\frac{D}{F}} \right)}} \cdot \left( {\frac{e}{{e_{2} }}} \right)^{{13\,\left( {\frac{D}{F}} \right)}} & = \frac{C}{2} \cdot F{\kern 1pt} \cdot \left( {a_{2} } \right)^{{ - {\kern 1pt} \frac{13}{2}}} \cdot \int {dt} , \Rightarrow e = e_{2} \cdot \left( {\frac{13C \cdot D}{2} \cdot \left( {a_{2} } \right)^{{ - \frac{13}{2}}} \cdot \Delta t} \right)^{{\,\frac{F}{13D}}} \\ \end{aligned}$$

(4.4)

where the scale-factor a ₂ should be given by the initial conditions according to the assumption e → 0.

Conclusion

We could obtain from the 1st of formulae (4.3) in “Appendix” above (as a first approximation), using the assumption e → 0 and then using the formulae (6.1):

$$\begin{aligned} \frac{da}{de}{\kern 1pt} = 2\frac{{\left( {D + E \cdot e^{2} } \right)}}{(1/a) \cdot F \cdot e}, \Rightarrow \frac{da}{a} \cong 2\left( {\frac{D}{F}} \right)\frac{de}{e}, \hfill \\ \Rightarrow \int {\left( {\frac{F(a)}{a \cdot D(a)}} \right)da} \cong \ln \left( {\frac{{e^{2} }}{{e^{2} (0)}}} \right), \Rightarrow \exp \left( {\int {\left( {\frac{F(a)}{a \cdot D(a)}} \right)da} } \right) \cong \left( {\frac{e}{e(0)}} \right)^{2} , \hfill \\ \end{aligned}$$

(6.2)

where

$$F(a) = \frac{{19k_{2}^{\kern 1pt} \,mR^{5} }}{{4\left( {2\left| {\sqrt {\frac{G(M + m)}{{a^{3} }}} -\upomega_{p} } \right|} \right)^{\upalpha} \cdot M}} - \frac{{7k_{{{\kern 1pt} 2}}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}},\quad D(a) = \frac{{k_{2} m \cdot R^{5} }}{{\left( {2\left| {\sqrt {\frac{G(M + m)}{{a^{3} }}} -\upomega_{p} } \right|} \right)^{\upalpha} \cdot M}}$$

for the definiteness, we will consider the case ω_p > n in the expressions above, so we should assume |ω_p − n| = ω_p − n. The left part of the Eq. (6.2) could be transformed as below ({a(0), e(0)} = const:

$$\begin{aligned} \exp \int {\left( {\frac{19}{4a} - \frac{{7k_{{{\kern 1pt} 2}}^{{{\kern 1pt} s}} M^{2} \cdot \left( {R_{s} } \right)^{5} }}{{k_{2} \cdot Q^{s} m^{2} R^{5} }}\frac{{\left( {2\left( {\upomega_{p} - \sqrt {\frac{G(M + m)}{{a^{3} }}} } \right)} \right)^{\upalpha} }}{a}} \right){\kern 1pt} da} \cong \left( {\frac{e}{e(0)}} \right)^{{{\kern 1pt} 2}} , \Rightarrow \hfill \\ \left( {\frac{a}{a(0)}} \right)^{{\frac{19}{4}}} \cdot \exp \left( { - \left( {\frac{{7k_{2}^{s} M^{2} \cdot \left( {R_{s} } \right)^{5} }}{{k_{2} \cdot Q^{s} m^{2} R^{5} }}\frac{{\left( {2\left( {\upomega_{p} - \sqrt {\frac{G(M + m)}{{a^{3} }}} } \right)} \right)^{\upalpha} }}{a}} \right)da} \right) \cong \left( {\frac{e}{e(0)}} \right)^{2} , \Rightarrow \hfill \\ \left( {\frac{a}{a(0)}} \right)^{{\frac{19}{4}}} {\kern 1pt} \cdot \exp \left( { - \frac{{7k_{2}^{s} M^{2} \cdot \left( {R_{s} } \right)^{5} }}{{k_{2} \cdot Q^{s} m^{2} R^{5} }}\int {\left( {\frac{{\left( {2\upomega_{p} } \right)^{\upalpha} \cdot \left( {1 -\upalpha\frac{{\sqrt {\frac{G(M + m)}{{a^{3} }}} }}{{\omega_{p} }} + \frac{{\upalpha \cdot (\upalpha - 1)}}{2}\left( {\frac{{\sqrt {\frac{G(M + m)}{{a^{3} }}} }}{{\upomega_{p} }}} \right)^{2} - ( \ldots )} \right)}}{a}} \right)} da} \right){\kern 1pt} \cong \left( {\frac{e}{e(0)}} \right)^{2} , \hfill \\ \end{aligned}$$

(6.3)

where we should restrict our approximation of the expression above for the range of the chosen parameter α = 0.16 ÷ 0.4 (in case of Saturn) by a proper way:

$$\begin{aligned} \left( {1 -\upalpha\frac{{\sqrt {\frac{G(M + m)}{{a^{3} }}} }}{{\upomega_{p} }} + \frac{{\upalpha \cdot (\upalpha - 1)}}{2}\left( {\frac{{\sqrt {\frac{G(M + m)}{{a^{3} }}} }}{{\upomega_{p} }}} \right)^{{{\kern 1pt} 2}} - ( \ldots )} \right) \hfill \\ \quad \cong 1 -\upalpha\frac{{\sqrt {\frac{G(M + m)}{{a^{3} }}} }}{{\upomega_{p} }} + \frac{{\upalpha \cdot (\upalpha - 1)}}{2}\left( {\frac{{\sqrt {\frac{G(M + m)}{{a^{3} }}} }}{{\upomega_{p} }}} \right)^{{{\kern 1pt} 2}} - \frac{{\upalpha \cdot (\upalpha - 1) \cdot (\upalpha - 2)}}{6}\left( {\frac{{\sqrt {\frac{G(M + m)}{{a^{3} }}} }}{{\upomega_{p} }}} \right)^{3} . \hfill \\ \end{aligned}$$

(6.4)

Indeed, an inequality below is valid for all the satellites in case of Saturn:

$$\left( {\frac{{\sqrt {\frac{G(M + m)}{{a^{3} }}} }}{{\upomega_{p} }}} \right) \le \sim\frac{2}{3},$$

so, we obtain for the 4-th term of the Taylor series (in the left part of (6.4)) that even at optimal conditions it should be less than

$$+ \frac{{\upalpha \cdot (\upalpha - 1) \cdot (\upalpha - 2) \cdot (\upalpha - 3)}}{24}\left( {\frac{2}{3}} \right)^{{{\kern 1pt} 4}} \cong + \frac{0.4 \cdot ( - 0.6) \cdot ( - 1.6) \cdot ( - 2.6)}{24}\left( {\frac{2}{3}} \right)^{4} \cong - 0.008$$

Using the Taylor decomposition of the expression (6.4) insofar, we obtain from (6.3) the appropriate invariant as below

$$\begin{aligned} \left( {\frac{a}{a(0)}} \right)^{{\frac{19}{4}}} \cdot \exp \left( {\begin{array}{*{20}c} {} \\ {} \\ \end{array} } \right. - \,\,\frac{{7\,k_{{{\kern 1pt} 2}}^{{{\kern 1pt} s}} \,M^{{{\kern 1pt} 2}} \, \cdot \,\left( {R_{\,s}^{\,} } \right)^{\,5} \,}}{{k_{{{\kern 1pt} 2}}^{\kern 1pt} \cdot Q^{\,s} \,m^{{{\kern 1pt} 2}} R_{\,}^{\,5} \,}} \hfill \\ \quad \cdot \int {\left( {\frac{{\left( {2\omega_{{{\kern 1pt} p}} } \right)^{{{\kern 1pt} \alpha }} \cdot \,\left( {1 - \alpha \frac{{\sqrt {\frac{G\,(M + m)}{{a^{{{\kern 1pt} 3}} }}} }}{{\omega_{{{\kern 1pt} p}} }} + \frac{\alpha \cdot (\alpha - 1)}{2}\left( {\frac{{\sqrt {\frac{G\,(M + m)}{{a^{{{\kern 1pt} 3}} }}} }}{{\omega_{{{\kern 1pt} p}} }}} \right)^{{{\kern 1pt} 2}} - \frac{\alpha \cdot (\alpha - 1) \cdot (\alpha - 2)}{6}\left( {\frac{{\sqrt {\frac{G\,(M + m)}{{a^{{{\kern 1pt} 3}} }}} }}{{\omega_{{{\kern 1pt} p}} }}} \right)^{{{\kern 1pt} 3}} } \right)}}{a}} \right)da} {\kern 1pt} \cong \left( {\frac{e}{e\,(0)}} \right)^{2} , \hfill \\ \quad \Rightarrow \left( {\frac{a}{a(0)}} \right)^{{\frac{19}{4}}} \cdot \exp \,\,\,\left( {\begin{array}{*{20}c} {} \\ {} \\ \end{array} } \right. - \,\,\frac{{7\,k_{{{\kern 1pt} 2}}^{{{\kern 1pt} s}} \,M^{{{\kern 1pt} 2}} \, \cdot \,\left( {R_{\,s}^{\,} } \right)^{\,5} \,}}{{k_{{{\kern 1pt} 2}}^{\kern 1pt} \cdot Q^{\,s} \,m^{{{\kern 1pt} 2}} R_{\,}^{\,5} \,}} \cdot \left( {2\omega_{{{\kern 1pt} p}} } \right)^{{{\kern 1pt} \alpha }} {\kern 1pt} \hfill \\ \quad \cdot \left( {\ln \,a + \alpha \frac{{\sqrt {G\,(M + m)} }}{{\omega_{{{\kern 1pt} p}} }}\frac{2}{3}\frac{1}{{a^{{\frac{3}{2}}} }} - \frac{\alpha \cdot (\alpha - 1)}{2}\frac{G\,(M + m)}{{(\omega_{{{\kern 1pt} p}} )^{{{\kern 1pt} 2}} }}\frac{1}{3}\frac{1}{{a^{3} }} + \frac{\alpha \cdot (\alpha - 1) \cdot (\alpha - 2)}{6}\left( {\frac{{\sqrt {G\,(M + m)} }}{{\omega_{{{\kern 1pt} p}} }}} \right)^{{{\kern 1pt} 3}} \frac{2}{9}\frac{1}{{a^{{\frac{9}{2}}} }} + const} \right) \cong \left( {\frac{e}{e\,(0)}} \right)^{2} , \hfill \\ \end{aligned}$$

where it is more than obvious that the terms below also should be neglected at the sufficiently large meaning of semi-major axis a (for the first approximation):

$$\begin{aligned} \alpha \frac{{\sqrt {G(M + m)} }}{{\omega_{p} }}\frac{2}{3}\frac{1}{{a^{{\frac{3}{2}}} }} - \frac{\alpha \cdot (\alpha - 1)}{2}\frac{G(M + m)}{{(\omega_{p} )^{2} }}\frac{1}{3}\frac{1}{{a^{3} }} + \frac{\alpha \cdot (\alpha - 1) \cdot (\alpha - 2)}{6}\left( {\frac{{\sqrt {G(M + m)} }}{{\omega_{p} }}} \right)^{3} \frac{2}{9}\frac{1}{{a^{{\frac{9}{2}}} }} \Rightarrow \hfill \\ 0.4\left( {\frac{2}{3}} \right)\left( {\frac{2}{3}} \right) - \frac{0.4 \cdot (0.4 - 1)}{2}\left( {\frac{2}{3}} \right)^{{{\kern 1pt} 2}} \frac{1}{3} + \frac{0.4 \cdot (0.4 - 1) \cdot (0.4 - 2)}{6}\left( {\frac{2}{3}} \right)^{3} \frac{2}{9}{\kern 1pt} \cong 0.18 + 0.018 + 0.004 < < \ln a \hfill \\ \end{aligned}$$

so, finally we obtain the proper invariant for the mutual dependence of the eccentricity e and the semi-major axis a as below (as the first approximation)

$$\left( {\frac{e}{e(0)}} \right) \cong \left( {\frac{a}{a(0)}} \right)^{H} , \Rightarrow a \cong a(0) \cdot \left( {\frac{e}{e(0)}} \right)^{{\frac{1}{H}}} ,H = \left( {\frac{19}{8} - \frac{{7k_{2}^{s} M^{2} \cdot \left( {R_{s} } \right)^{5} }}{{2k_{2} \cdot Q^{s} m^{2} R^{5} }} \cdot \left( {2\upomega_{p} } \right)^{\upalpha} } \right)$$

(6.5)

Using (6.5), we could obtain from the 2nd of Eqs. (4.1) and expressions (4.2), (6.2):

$$\frac{de}{dt} = \frac{C}{{2(a(0))^{{\frac{13}{2}}} \cdot \left( {\frac{e}{e(0)}} \right)^{{{\kern 1pt} \frac{13}{2H}}} }} \cdot F \cdot e,$$

$$\begin{aligned} C{\kern 1pt} = \pm 3\sqrt {G(M + m)} ,F = \frac{{19k_{2} mR^{5} }}{{4\left( {2\left( {\upomega_{p} - \sqrt {\frac{G(M + m)}{{(a(0))^{3} \cdot \left( {\frac{e}{e(0)}} \right)^{{\frac{3}{H}}} }}} } \right)} \right)^{\upalpha} \cdot M}} - \frac{{7k_{2}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} \,m\,}}, \Rightarrow \hfill \\ \int {\frac{{\left( {\upomega_{p} - \sqrt {\frac{G(M + m)}{{(a(0))^{3} \cdot \left( {\frac{e}{e(0)}} \right)^{{\frac{3}{H}}} }}} } \right)^{\alpha } \cdot \left( {\frac{e}{e(0)}} \right)^{{{\kern 1pt} \left( {\frac{13}{2H} - 1} \right)}} }}{{\left( {\frac{{19k_{2} mR^{5} }}{{4M \cdot 2^{\upalpha} }} - \frac{{7k_{2}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{Q^{s} m}} \cdot \left( {\upomega_{p} - \sqrt {\frac{G(M + m)}{{(a(0))^{3} \cdot \left( {\frac{e}{e(0)}} \right)^{{\frac{3}{H}}} }}} } \right)^{\upalpha} } \right)}}d\left( {\frac{e}{e(0)}} \right)} = \left( {\frac{C}{{2(a(0))^{{\frac{13}{2}}} }}} \right) \cdot \int {dt} \hfill \\ \end{aligned}$$

(6.6)

where (6.6) yields (as the first approximation by Tailor series):

$$\int {\frac{{\left( {\frac{e}{e(0)}} \right)^{{\left( {\frac{13}{2\,H} - 1} \right)}} }}{{\left( {\frac{{19k_{2} mR^{5} }}{{8\left( {2\upomega_{p} } \right)^{\upalpha} \cdot M}} \cdot \left( {1 +\upalpha\sqrt {\frac{G(M + m)}{{(a(0))^{3} \cdot (\upomega_{p} )^{2} }}} \cdot \left( {\frac{e}{e(0)}} \right)^{{ - \frac{3}{2H}}} } \right) - \frac{{7k_{2}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{2Q^{s} m}}} \right)}}d\left( {\frac{e}{e(0)}} \right)} \cong \left( {\frac{C}{{(a(0))^{{\frac{13}{2}}} }}} \right) \cdot \int {dt}$$

(6.7)

The left part of Eq. (6.7) could be simplified by a proper change of variables:

$$\begin{aligned} u = \left( {\frac{e}{e(0)}} \right)^{{ - \frac{3}{2H}}} \Rightarrow \left( {\frac{e}{e(0)}} \right) = u^{{ - {\kern 1pt} \frac{2H}{3}}} \Rightarrow d\left( {\frac{e}{e(0)}} \right) = - \left( {\frac{2H}{3}} \right) \cdot u^{{\left( { - \frac{2H}{3} - 1} \right)}} du \Rightarrow \hfill \\ \int {\left( {\frac{{u^{{\frac{2H}{3}\left( {1 - \frac{13}{2H}} \right)}} }}{{\left( {\frac{{19k_{2} mR^{5} }}{{8\left( {2\upomega_{p} } \right)^{\alpha } \cdot M}} \cdot \left( {1 +\upalpha\sqrt {\frac{G(M + m)}{{(a(0))^{3} \cdot (\upomega_{p} )^{2} }}} \cdot u} \right) - \frac{{7k_{2}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{2Q^{s} m}}} \right)}}} \right) \cdot u^{{\left( { - \frac{2H}{3} - 1} \right)}} du} = - \left( {\frac{3}{2H}} \right) \cdot \left( {\frac{C}{{(a(0))^{{\frac{13}{2}}} }}} \right) \cdot \int {dt} , \hfill \\ \int {\left( {\frac{{u^{{\;{\kern 1pt} \left( {\; - \;\frac{16}{3}} \right)}} }}{{\left( {\frac{{19k_{2} mR^{5} }}{{8\left( {2\upomega_{p} } \right)^{\alpha } \cdot M}} \cdot \left( {1 +\upalpha\sqrt {\frac{G(M + m)}{{(a(0))^{3} \cdot (\upomega_{p} )^{2} }}} \cdot u} \right) - \frac{{7k_{2}^{s} M \cdot \left( {R_{s} } \right)^{5} }}{{2Q^{s} m}}} \right)}}} \right)du} = - \left( {\frac{3}{2H}} \right) \cdot \left( {\frac{C}{{(a(0))^{{\frac{13}{2}}} }}} \right) \cdot \int {dt} , \hfill \\ \end{aligned}$$

(6.8)

where the left part of the last of Eqs. (6.8) could be presented as below

$$\begin{aligned} \int {\frac{{u^{{\left({- \frac{16}{3}} \right)}}}}{{\left({\Lambda \cdot \left({1 \pm \left({\frac{\upalpha C}{{3\upomega_{p} \cdot (a(0))^{{\frac{3}{2}}}}}} \right) \cdot u} \right) - \varXi} \right)}}} du = - \left({\frac{3}{{2H \cdot (a(0))^{{\frac{13}{2}}}}}} \right) \cdot \int {dt}, \hfill \\ \Lambda = \frac{{19k_{2} mR^{5}}}{{8\left({2\upomega_{p}} \right)^{\upalpha} \cdot M}},\quad \varXi = \frac{{7k_{2}^{s} M \cdot \left({R_{s}} \right)^{5}}}{{2Q^{s} m}} \hfill \\ \end{aligned}$$

(6.9)