Climate Vibrations and Chaotic Behavior of Earth’s Inclination in the Laskar’s Model of a Moonless Earth

The model of the moonless Earth, introduced by J. Laskar, has the form of a non-autonomous Hamiltonian system of differential equations for two variables: the cosine of the angle of inclination and the longitude of the axis of rotation of the Earth. The system describes the rotational dynamics of the Earth under the influence of the sun and planets. Earth perturbations from other planets of the solar system are considered periodic and are taken into account using the first four terms of the Fourier expansion of the corresponding part of the Hamilton function with known amplitudes and frequencies. The initial inclination of the Earth is considered as a parameter of the problem. The system was numerically integrated over a time period of 18 million years for various values of the initial inclination from 0 to 180 degrees. Three chaotic gaps of the initial inclination were found. During the bifurcation study, singular points were found and special segments of the non-autonomous system were obtained. A bifurcation diagram of the system is constructed by the initial inclination parameter. Poincare cartographic maps are constructed. The system is written in variations on the initial conditions for the Laskar system, and with its help the dependences of the problem parameter of the senior Lyapunov exponent and the averaged MEGNO indicator are calculated. The results confirm the presence of three chaotic and one regular region of variation of the bifurcation parameter of the problem.


Introduction
In recent years the dynamics of the axis of rotation of the moonless Earth is an important task for assessing potential climate changes. We can apply a model of a moonless Earth for many planetary systems for which the masses and orbit parameters are well known. But, unfortunately, the results will not be sufficiently informative due to the lack of data on the observation of the physical properties of individual planets, such as rotational speed and dynamic ellipticity. Consequently, detailed studies should expect significant breakthroughs in the observational technique in astronomy.
However, the implications of modern research for extrasolar planetary systems are clear: the lack of a large mass ratio of the moon and the earth should not be seen as evidence of extreme climate change. Thus the characteristic chaotic diffusion rate can be slow enough not to limit the long-term suitability for life on the planet.

Formulation of the problem
We will study the evolution of the axis of rotation of a moonless Earth in the context of a simplified perturbed system.
The Hamiltonian for the model of the moonless Earth has the form: where -is the longitude of the axis of rotation, = cos( ), -is the slope, а , -are an approximately constant parameters.
, и -are the coefficients presented in Table 1.

Finding of the singular points
Find the singular points of the system (3) from conditions ̇= 0, ̇= 0: Since the system is not autonomous, we look for singular points ( , ) first for fixed .
So we found a lot of pairs ( , ). Examples for pairs ( 1 , 1 ), ( 2 , 2 ), ( 3 , 3 ) and ( 4 , 4 ) can be seen in the figures 1 -4.    To determine the type of the singular point, we compose the Jacobi matrix for our case and call it matrix A. Next, we find the eigenvalues of matrix A at our point ( , ). To do this, we calculate the determinant of the matrix А1 = А − * , where 1,2 eigenvalues, Е -is the unit matrix of dimension 2х2. Thus, we obtain (in our case) a quadratic equation whose roots will be our eigenvalues. The found eigenvalues have purely imaginary roots of different signs.

Bifurcation diagram
Now let's move on to the bifurcation diagram. To find it, we use the system 2 or 3 in expanded form. We take the initial values from the work of Laskar: The remaining coefficients are taken from the table 1 and equations 2 and 3. We will calculate for time = 18 ⋅ 10 6 . We will calculate for time: where 0 -initial time, а n -split. The result of solving a system of differential equations with initial conditions can be seen in the figures 5-8. Moreover, for different and charts will vary.

Lyapunov exponent
The next step in our study will be finding the Lyapunov exponent.
Consider the variations of the system solutions (2.5) with initial conditions: Then the system will take the form:

#(10)
Lyapunov exponent is calculated by the formula: The initial conditions will be the following equalities: The necessary coefficients are taken from the table. 1, equalities 2 и 3 will be added to them. = 18 * 10 6 , а ε = (20 * )/180. The result of the first run can be seen in the figure. 9.

Mean Exponential Growth of Nearby Orbit
Along with the general Lyapunov variational method, a method called MEGNO or Mean Exponential Growth of Nearby Orbit is widely used. Its essence is to calculate the average parameter ( ). This method is very popular in the scientific community, because it can be used to separate chaotic and regular modes of movement.
Parameter ( ) is a time-weighted integral form: where δ -variational vectors.
We write the average value in the form: Due to the fact that we have the task of numerical modeling, it would be more logical to replace the integral relations with the corresponding differential relations. Then we get: Then the parameters ( ) and � ( ) will take the form: Finally, the system for MEGNO will look like: As initial conditions, we take the previously mentioned equalities 12 and 13. Equals 25-28 will be added to them:

Integration of the Laskar system
Now back to the initial system 3. Solving it, we get the following graphs with dependencies:

Poincare mappings
Now we turn to the Poincaré map. To do this, make the transition ( , ) → (ε, ).
Having the initial conditions = cos (ε) find : ε = arccos ( ) and considering �1 − 2 = sin( ), we get system: Substitute the values (not forgetting the initial conditions) and rewrite the system completely:

#(30)
In figures 17-20 the Poincare map is shown at different approximations. Separately, it is worth considering that as a period t coefficient 3 adopted, as the maximum parameter.

Conclusion
In this work, we examined the Laskar model of the moonless Earth. Special segments of the nonautonomous Laskar system with singular points of the "center" type were found. The evolutions of the inclination and longitude of the axis of rotation over the interval 10 7 years are obtained.
Three chaotic regions of the evolution of inclination on the axis of the initial inclination parameter are found. The bifurcation diagrams of the system of the system are constructed for various parameters , 0T. Plots of the Lyapunov exponent and the averaged MEGNO exponent are plotted depending on the initial inclination. A stroboscopic Poincare mapping corresponding to a chaotic attractor is constructed. The presence of a chaotic attractor explains the absence of catastrophic climate changes in the chaotic region of the evolution of inclination.