New Long-Term Climate Oscillations ()

Joseph J. Smulsky^{}

Institute of Earth’s Cryosphere, Tyum. SC of SB RAS, Federal Research Center, Tyumen, Russia.

**DOI: **10.4236/ns.2021.138028
PDF HTML XML
112
Downloads
737
Views
Citations

Institute of Earth’s Cryosphere, Tyum. SC of SB RAS, Federal Research Center, Tyumen, Russia.

The astronomical theory of climate change is based on the solution of differential equations describing Earth’s orbital and rotational motions. The equations are used to calculate the change in insolation over the Earth’s surface. As a result of the author’s solution of the orbital problem, the periods and amplitudes of Earth-orbit variations and their evolution have been refined. Unlike previous studies, the equations of Earth’s rotational motion are solved completely. The Earth’s rotational axis precesses relative to a direction different from the direction of the orbit’s axial precession, and oscillates with periods of half a month, half a year and 18.6 years. Also, its oscillations occur with irregular periods of several tens of thousands of years and more. All these motions lead to oscillations of the obliquity in the range of 14.7° to 32.1°, which prove to be 7 - 8 times larger than obtained by a previous theory. In the same proportion, the Earth’s insolation oscillations increase in amplitude, with insolation extremes occurring in other epochs than those in the previous theory. The amplitudes and the onset times of the extremes correlate with known paleoclimate changes. Thirteen insolation periods of paleoclimate variation over an interval of 200 thousand years are identified. From the insolation evolution calculated over a time interval of 1 million years, 6 climate gradations from very cold to very warm are identified.

Share and Cite:

Smulsky, J.J. (2021) New Long-Term Climate Oscillations. *Natural Science*, **13**, 354-371. doi: 10.4236/ns.2021.138028.

1. INTRODUCTION

Long-term climate oscillations are analyzed in the Astronomical theory of climate change (or alternately named, Astronomical theory of ice ages). The theory is based on the solution of the following three problems: 1) what are the changes in the Earth’s orbit? 2) what are the changes in the Earth’s axis of rotation? 3) what are the changes in the amount of solar radiation over the Earth’s latitude based on the first two changes? The original version of the theory was developed by M. Milankovitch [1] in the first quarter of the 20^{th} century. Subsequently, the proposed approach was improved by other researchers [2 - 6]. At the end of the 20^{th} century, activities aimed at revisiting the above problems were initiated [7]. As the result of a more precise solution to the problem of Earth’s rotational motion, a second version of the climate change theory succeeded in explaining the factors that cause long-term climate oscillations.

This paper presents the main results of the new Astronomical theory of climate change, how they were obtained and their reliability.

2. METHOD

This paper presents the results of decades of theoretical studies. The foundations of mechanics were analyzed, unsubstantiated hypotheses were eliminated and mechanical operations were formulated with regard to their actual content [8]. The differential equations of the orbital motion of planets in the Solar system were redefined, a new method for solving them was developed, and the Galactica program was created for their numerical solution [9 - 11]. Differential equations for orbital and rotational motion were derived in a new way and a method for analyzing, solving and integrating them was developed. This evolved into a program for their numerical solution [11 - 13]. Orbital motion was considered in the coordinates associated with the fixed equator of Earth, and rotational motion was investigated in coordinates associated with a fixed plane of the Earth’s orbit. Transformations of the parameters of orbital motion relative to the coordinates associated with the plane of the moving equator are derived [14]. The equations for the irradiation of the Earth by the Sun, depending on the parameters of its orbital and rotational motions, are derived in a new way. A program for calculating the Earth insolation has been developed [15]. The numerical integration of these equations on supercomputers has been performed. These tasks are complex and laborious. For example, the integration of the equations of orbital motion over 100 million years took two years of machine time. The analysis of the solutions of these equations was carried out over several years.

These studies are unique, so all their stages were tested. Algorithms for checking the accuracy of solutions were introduced into the programs for integrating equations. For example, there are about two dozen in the Galactica program. All stages of the solution were checked, compared with the observations and conclusions of other researchers. Checking took longer than the solutions themselves. For example, verification of solutions to a problem on the rotation of the Earth was carried out over a period of three years.

3. RESULTS

3.1. Earth’s Motions and Their Variations

The Earth moves in an elliptical orbit around the Sun, which is located at one focus of an ellipse (Figure 1). The shortest Earth-Sun distance in perihelion is denoted by *R _{p}*, and the largest distance in aphelion by

With respect to motionless space, the Earth rotates around its axis,
$\stackrel{\to}{N}$ at an angular velocity of *ω _{E}* = 7.292115∙10

Figure 1. The Earth’s position in its orbit in 2025 during the days of vernal equinox (20.03), summer solstice (21.06), autumnal equinox (22.09), and winter solstice (21.12), and the time expressed in Earth motion during the days in spring (92.7 d), in summer (93.7 d), in autumn (89.9 d), and in winter (89.0 d): $\stackrel{\to}{N}$ is the Earth’s axis of rotation; ${\stackrel{\to}{M}}_{2}$ is a vector relative to which the axis $\stackrel{\to}{N}$ precesses in a period of 25.74 thousand years; $\stackrel{\to}{S}$ is the Earth’s rotational axis, and $\stackrel{\to}{M}$ is a vector relative to which the $\stackrel{\to}{S}$ axis precesses over a period of 68.7 thousand years [14 , 15 , 18 , 19].

is at maximum inclination to the Earth-Sun line. This is the reason the southern hemisphere is at maximum illumination at that time, and polar night approaches at high latitudes in the northern hemisphere. Since the time spent on reaching and leaving the extreme angles lasts for several days, these points are called respectively, the summer solstice day (21.06), and the winter solstice day (21.12).

The inclination of Earth’s axis, $\stackrel{\to}{N}$ to the orbital axis, $\stackrel{\to}{S}$ leads to the variation in duration of sunlight, both during the year and on the same day at different latitudes. On summer solstice day, (Figure 1, 21.06), we have a polar day in the whole region between the North Pole and the Arctic Circle. Then, as the latitude decreases, the day gets shorter, reaching a 12-hour duration at the equator, and we have a polar night established below the Antarctic Circle. Contrary to this, on the day of winter solstice, (21.12), in the territory between the North Pole and the Arctic Circle we have a polar night. Then the day starts increasing in duration. At the equator, the day lasts for 12 hours, and a polar day sets in below the Antarctic Circle. As we approach the equinoctial points of 20.03, and 22.09, the difference between the days in latitude decreases in value. The day’s duration along all latitudes becomes identical at 12 hours.

As the Earth moves along its orbit, the alteration of seasons occurs. The duration of the seasons is defined by the Earth’s motion over certain orbital segments. From the vernal equinox day, 20.03, until the summer solstice day, 21.06, the duration of spring is 92.7 days. Over the *summer* segment, the duration is 93.7 days. Over the *autumnal* segment, the duration is 89.9 days. Over the *winter* segment, the duration is 89.0 days.

The Earth’s orbital and rotational motions define the variation in climate in the current epoch. However, the motions vary in time, and the climate undergoes change. The position of Earth’s orbit precesses in space. The Earth’s orbital axis
$\stackrel{\to}{S}$ (Figure 1), rotates, or in other words, it precesses about the direction of
$\stackrel{\to}{M}$ , which is motionless in space. The precession proceeds clockwise with a period of 68.7 thousand years. In a clockwise direction, the Earth’s axis
$\stackrel{\to}{N}$ precesses about the direction of
${\stackrel{\to}{M}}_{2}$ , also motionless in space. The precession period is 25.74 thousand years. Besides this, the axes
$\stackrel{\to}{S}$ and
$\stackrel{\to}{N}$ execute oscillations, each with respect to its own precessional axis
$\stackrel{\to}{M}$ and
${\stackrel{\to}{M}}_{2}$ , respectively. In addition, the shape of the orbit (its eccentricity, whose value varies from 0 to 0.064; the current value being equal to 0.016) and the perihelion position, both undergo variations. Today, the perihelion is over the *winter* segment (Figure 1), when winter sets in the northern hemisphere. Since the Earth’s orbital perihelion rotates in anticlockwise direction at a mean period of 147 thousand years, its position in other epochs can be at any point in the Earth’s orbit.

The changes of the Earth’s orbital and rotational motions lead to changes in its climate. Below, the latter changes are considered in more detail.

3.2. Evolution of the Earth’s Orbital Motion

The evolution of the Earth’s orbital motion is considered in a motionless frame, *xyz* whose origin is located at the center, *O* of celestial sphere 1 (Figure 2). Note that, depending on the particular problem of interest, the point *O* can be located either at the center of mass of the Solar system, at the center of the Sun, or at the center of the Earth. As a result of the interaction of Solar-system bodies, the Earth’s equatorial plane 2 and its orbital plane 3 both alter their positions, which are denoted with digits 2 and 3 in the epoch *Т*_{0}, and with digits 4 and 5 in the epoch* Т,* respectively. At epoch* Т*_{0}, we consider the epoch of the year 2000.0.

Since the annual motion of the Sun over celestial sphere 1 with respect to the Earth, proceeds along circles 5 or 3, the planes of the circles are additionally called, the planes of moving and motionless ecliptic, respectively. The frame *xyz* is related to the plane of the motionless Equator 2. The moving plane 5 of the Earth’s orbit is defined by the angle *φ*_{Ω} = *γ*_{0}*γ*_{2} of the ascending-node position* γ*_{2} and by the inclination of the plane *i*.

The Earth moves around the Sun along an open trajectory which is close to the shape of an ellipse. At one point in the trajectory, the perihelion, the Earth approaches the Sun to the shortest distance *R _{p}*, and at the opposite point, the aphelion, it moves away from the Sun to the largest distance

The oscillations of the angles *φ*_{Ω} and *I* reflect the rotation process of the orbit’s axis
$\stackrel{\to}{S}$ (Figure 2), that is, the normal to the orbit plane, with a period of 68.7 thousand years about the motionless vector
$\stackrel{\to}{M}$ [20]. The latter vector is the sum of the angular momentum vectors of all bodies in the Solar system. The plane normal to
$\stackrel{\to}{M}$ is called the Laplace plane. The angles that define the position of the Laplace plane with respect to the motionless equatorial plane *xOy* in radians are equal to *i _{M}* = 0.401834 and

Figure 2. Parameters of the Earth’s orbit and axis in the motionless equatorial *xyz* and motionless ecliptic *x _{e}y_{e}z_{e}* coordinate systems. 1—the celestial sphere; the planes in the epoch of

in clockwise direction. In addition to the precession, the orbit’s axis
$\stackrel{\to}{S}$ executes oscillations at periods of 97.35 thousand years, 1.164 million years, and 2.32 million years, relative to the vector
$\stackrel{\to}{M}$ . During those oscillations, the angle between the vectors
$\stackrel{\to}{S}$ and
$\stackrel{\to}{M}$ never exceed a value of 2.94˚. All those motions are manifested in the behavior of angles *φ*_{Ω} and *i* in Figure 3.

Thus, the evolution of the Earth’s orbit proceeds as a result of the following four motions: 1) precession of the orbit’s axis
$\stackrel{\to}{S}$ ; 2) oscillations of the orbit’s axis
$\stackrel{\to}{S}$ ; 3) oscillations of the orbit’s eccentricity*e*; and 4) rotation of the orbit in its own plane (perihelion rotation).

Figure 3 shows the evolution of the Earth’s orbital parameters *e*, *φ*_{Ω}, *i* and *φ _{p}* over a span of one million years. The shortest eccentricity oscillation period is

3.3. Evolution of the Earth’s Rotational Motion

Evolution of the Earth’s rotational motion is treated in the motionless frame *x _{e}y_{e}z_{e}* (Figure 2) connected

Figure 3. Evolution of Earth’s orbital parameters over a time interval of 1 million years: *е*—eccentricity; *φ*_{Ω}—angle of the orbit’s ascending node; *i*—orbit’s inclination; *φ _{р}*—perihelion angle;

with the Earth orbit’s motionless plane, 3. The inclination angle and the precession angle *ψ *= *γ*_{0}*γ*_{1} both define the position of equator, 4 moving with respect to the orbit’s motionless plane, 3. The precession angle *ψ* decreases toward future in an oscillatory manner at a mean rate of
${\stackrel{\dot{}}{\psi}}_{m}=-\text{2}\cdot \pi /{P}_{pr}$ according to a linear law. Here *P _{pr}* = 25.738 thousand years is the period of the Earth axisial precession averaged over a time interval of 1 million years [7 , 13 , 14]. The negative sign of
${\stackrel{\dot{}}{\psi}}_{m}$ implies that the unit vector
$\stackrel{\to}{N}$ of the Earth’s rotation axis precesses in a clockwise direction (Figure 2). Precession of the Earth’s rotational axis,
$\stackrel{\to}{N}$ proceeds around the second motionless vector
${\stackrel{\to}{M}}_{2}$ ; the angles defining the orientation of the plane normal to this vector with respect to the plane

The variation of the difference,
$\Delta \psi =\psi -{\psi}_{a}$ for the precession angle, *ψ* (here,
${\psi}_{a}={\psi}_{0}+{\stackrel{\dot{}}{\psi}}_{m}\cdot T$ is the change of the precession angle in the process proceeding at the mean rate
${\stackrel{\dot{}}{\psi}}_{m}$ ) over a time interval of 1 million years, is shown in Figure 4. Over the latter time interval, the quantity, Δ*ψ* varies from −0.184 to 0.233 radian, so that the full oscillation swing amounts to 0.417 radian.

The inclination difference,
$\Delta \theta =\theta -{\theta}_{0}$ in Figure 4, is given with respect to the initial angle *θ*_{0} = 0.40904645. The quantity *Δθ* oscillates similar to, Δ*ψ**,* yet in a narrower range from −0.0845 to 0.0855, so that the oscillation range here amounts to 0.17 rad. Thus, the oscillation amplitude for angle *θ* is 2.45 times smaller than that for angle *ψ*. Besides, the oscillations of Δ*θ* do not coincide in phase with the oscillations of Δ*ψ*; they are shifted along the time axis *T* by −7.5 thousand years.

3.4. Evolution of the Earth’s Orbital Motion Relative to Its Rotational Motion

The orbital-motion parameters *i*, *φ _{Ω}* and

Figure 5 shows the variation of *ε* over five different time intervals *In* [13]. Over short time intervals, the oscillations of *θ* and *ε* are roughly identical. Indicated in the graphs are the main oscillation periods *T _{ni}* and amplitudes (

Over the time interval *In* = 0.1 year, half-monthly oscillations are observed and diurnal oscillations with period *T _{n}*

Over the time interval *In* = 100 years, it is seen that the calculated obliquity *ε* 1) oscillates about its mean value, 2) received by S. Newcomb [16] and J. Simon *et al.* [17]. The oscillation amplitude *ε _{a}*

Figure 4. Evolution of the Earth’s rotational motion over a time interval of 1 million years. The differences; the precession Δ*ψ* and inclination Δ*θ* angles are given in radians.

Figure 5. The dynamic of the obliquity *ε* (in radians) over five time intervals *I**n*: *yr*—year; Δ*θ* ≈ *ε *− *ε*_{0}; *ε*_{0}—the obliquity at the initial epoch of December 30, 1949. *T _{n}*

3.5. Evolution of the Earth’s Obliquity and Insolation over a Span of 1 Million Years

As it is evident from Figure 5, over the time interval of *In* = 10 thousand years, a coincidence of the new obliquity *ε* 1 with the data 2 and 3 yielded by the first version of the Astronomical theory of climate change [1 - 6] is observed over a span of 2000 years.

It should be noted that these authors solved the problem in question over a large time interval: Sharaf and Budnikova [3] for 30 million years, and Laskar *et al.* [6] for an even longer period. We showed [18] that their results do not fundamentally differ from the results of Milankovitch [1], Berger and Loutre [4], Edvardsson *et al.* [5]. Therefore, we compare our results with these authors who are typical representatives of the former Astronomical theory.

As can be seen from Figure 5, after 2000 years, obliquity, *ε* calculated within the new version 1 of the theory shows clear deviation. As it is seen from the graphs of Figure 6, over the time interval of 1 million years the oscillations of *ε* as yielded by the second version of the theory proceed in the range of 14.7˚ to 32.1˚, whereas the same range in the previous theory was from 22.08˚ to 24.45˚; in other words, the range of oscillations in the second version of the theory proves to be seven times greater.

This difference is due to the fact that in the second version of Astronomical theory, the Earth’s rotation problem was treated in full, without simplifications. The solution of this problem and various checks of obtained data were analyzed at length in the publications [7 , 13 , 21] and will be covered in the subsection 4.3.

The amount of solar radiation reaching the Earth’s surface, also called the Earth’s insolation, is defined by the parameters *e*, *ε* and *φ _{pγ}*. Figure 6 gives a comparison of the insolation changes
${Q}_{s}^{\text{65}N}$ occurring during the summer caloric half-year at the 65-deg northern latitude in the second version of the theory, (curve 1) [21] with the changes as calculated by the first version of the theory [6]. Here, the amplitude of insolation oscillations is also seven times greater than that in the previous theory. Besides, the insolation extremes occur at other times, and the oscillation periods are different. Note that the astronomical

Figure 6. Evolution of the obliquity *ε* and that of summer half-year insolations
${Q}_{s}^{\text{65}N}$ and *I* over a time interval of 1 million years. Comparison of results yielded by the second version of the Astronomical theory of climate change (curve 1), with the results of the first version of this theory (curve 2), demonstrated using, as an example, the work by J. Laskar *et al.* [6].*
${Q}_{s}^{\text{65}N}$ —insolation in GJ/m ^{2} over the summer caloric half-year at 65-deg northern latitude; I—insolation at the equivalent latitudes over the summer caloric half-year at 65-deg northern latitude.* Indicated in degrees are the maximum and minimum values of

summer and winter half-years measured from the vernal equinox day to the autumnal equinox day and visa versa differ in duration for different epochs. That is why it is caloric half-year, equal in duration, that are considered here.

In order to compare climates in other epochs with the current climate, we consider the insolations at equivalent latitudes *I*. For calculating of *I*, we consider the Earth’s latitude *φ* characterized by receiving same amount of summer solar radiation, *Q _{s}* as in the current epoch. Figure 6 shows the insolation oscillations at equivalent latitudes

3.6. Variation of Insolation over the Earth’s Latitude

Figure 7 shows the variation over the latitude *φ* of the year *Q _{T}*, and during caloric half-years summer

The summer insolation *Q _{s}* (Figure 7, dashed lines) in contemporary epoch, 1 exhibits a minimum value on the poles, and it reaches a maximum value at the tropics

The winter insolation *Q _{w}* (Figure 7) on the poles is zero, and it monotonically increases in the equatorial region. In the equatorial region, the insolation

The annular insolation *Q _{T}* (Figure 7) monotonically increases from the poles toward the equator. At the equator, the annular insolation exhibits a maximum, with the annular insolation being symmetrical with respect to the equator. In other words, the amounts of heat per year are identical in both hemispheres. From cold epoch 3

Figure 7. The distribution over latitude *φ* of insolations over three epochs. (*Q _{s}*, is summer half-year insolation;

and at the latitude, *φ* = 45˚ the annular insolation experiences no changes. In the equatorial region, the changes of *Q _{T}* are reciprocal to its changes at the high latitudes: in cold epoch 3, the amount of heat per year exceed that in the warm epoch. In the latter situation, the change of insolation

3.7. Periods and Gradations of Earth’s Climate Changes

Over the previous interval of 200 thousand years (see Figure 8), 13 climatic periods *O _{I}*, 1

Also, the following gradations of the warm and cold climate were introduced (Figure 8): moderately warm, warm, and extremely warm climate levels, and moderately cold, cold, and extremely cold climate levels. During the past period of 1 million years (see Figure 9), the Earth has experienced six extremely cold (e.c.) periods and four extremely warm (e.w.) periods. The total number of cold (c.) and warm (w.) periods was 16 warm periods and 16 cold periods. Other periods were moderately cold (m.c.) ones and

Figure 8. Insolation periods *O _{I}*, 1

Figure 9. Climate levels over a time interval of 1 million years: 1—mean insolation *Q _{sm}*; 1

moderately warm (m.w.) ones. Besides, there were nine moderate-climate (m.cl.) periods, which included both cooling and warming phases.

4. DISCUSSION

4.1. Statements of the Problems and the Differences among Them

In the previous Astronomical theory of climate change, variations of Earth’s orbital elements were received on the basis of the theory of secular perturbations. This is an approximate analytical method for solving the interaction problem for Solar-system bodies. In that theory, the changes of the equator plane 4 (Figure 2) were analyzed approximately. With respect to this plane, the angle *ε* of the inclination of the orbital plane 5 and the perihelion position (being, in our notation, the angle *φ _{pγ}* in Figure 2) were determined.

In the new Astronomical theory the interaction problem for Solar-system bodies (for simplicity, we will call this the orbital problem) was solved, using no simplifications, with a high-precision numerical method using a specially developed Galactica system [9 - 11]. In the latter case, we consider a change in the orbital plane (item 5 in Figure 2) relative to the fixed space represented by the equatorial plane 2 on a certain epoch. Here, the inclination, *i* of the orbital plane 5 and the ascending-node angle *φ*_{Ω} differ from the angles figuring in the theory of secular perturbations. As it was noted above, in this theory the angle *ε* between the moving equatorial plane 4 and the moving orbital plane 5 was used.

Unlike in the theory of secular perturbations, we also numerically solved a second problem; the one about the Earth’s rotation, governed by its own differential equations. As a result, we obtain the laws of variation of the inclination angle *θ* and the precession angle *ψ* of the moving equatorial plane 4 with respect to the motionless orbital plane 3 (Figure 2).

A third complex geometric problem was then solved analytically for determining both the obliquity *ε* of the moving equatorial plane 4 with respect to the moving orbital plane 5 and the perihelion angle *φ _{pγ}*. The solution of the problem over short time intervals of the order of one thousand years is available in celestial mechanics. Over time intervals of several million years, during which those two planes and the orbit perihelion executed many irregular rotations in both directions, no solution was known.

A fourth problem, which yields the variation of the insolation versus the change of Earth’s orbital and its axis parameters, thus offering an insolation theory, was presented in its complete form by M. Milankovitch in the first quarter of the twentieth century. We have also solved this problem in a new way [15]. Here, a new mathematical algorithm for elliptic motion, more understandable to non-specialists and useful for computer calculations, was employed.

All equations, including the differential equations for the orbital and rotational motions, were derived using a novel approach. Since the resulting data were other quantities in different representations, new methods for their analysis were developed. All that activity was accompanied with the development of computer codes written in various programming languages.

4.2. Adequacy of the Solution of the Orbital Problem

In connection with the new solutions of all four problems, at due stages, all the problems were subjected to checking. For checking the adequacy in the solution of the orbital problem, nine trustworthiness criteria were developed. Some of those criteria were included into the Galactica program; that is why control was exercised right in the course of solving the problem. For all bodies having an observation base (these are the planets from Mercury to Neptune, and also the Moon, the solutions over the time interval of several thousand years were compared with the secular variations of orbital parameters. The coincidence was found to be excellent [20 , 23].

Over intervals of hundred thousand and million years, the orbital parameters were compared to the results obtained by the previous researchers [3 , 4 , 6]. The data were also found to be coincidental. Each of subsequent authors took into account the experience gained by previous researchers and made the theory of secular perturbations more refined. The theory was confirmed by making necessary comparisons. The later the works were published, the longer was the time interval over which the corresponding results were coincident with our data [20].

As it was noted previously, the perturbation theory is an approximate method for solving the orbital problem. After 20 million years, the solution obviously started departing from the actual data: the orbits of individual planets started increasing in size and, later, the planets could leave the Solar system [24]. We have solved the orbital problem over a time interval of 100 million years. All the orbital parameters of the planets and Moon executed steady oscillations, and there existed no tendency toward the change of those oscillations [20].

4.3. Adequacy of the Solution of the Earth’s Rotation Problem

Points concerning the reliability of the solution of the Earth’s rotation problem were analyzed in detail in publications [7 , 13 , 14]. Within the adopted solution method, all necessary checks were performed. For instance, the problem was solved in succession, implying the action of one of ten bodies (the Sun, the eight planets, and the Moon) [12]. The obtained oscillation periods of the Earth’s axis were confirmed with general conclusions made on the basis of the theorem of change in angular momentum and also, with the results of other authors [25]. With the actions due to all bodies, the problem was solved over different time intervals; as it was shown previously in sec. 4 and sec. 5, the obtained data proved to be coincident with observations. Integration of the equations over a time interval of 200 thousand years was performed with different initial conditions and integration steps. The latter resulted in no changes of the oscillation periods, oscillation amplitudes, and the onset times of the extremes.

From the graph with the interval *I _{n}* = 100 yr in Figure 5, one can see that the middle of the oscillations of the calculated 1 obliquity

Following 2.5 thousand years, solutions obtained by previous researchers also started departing from the linear trend. From this time on, differences between solutions 2 and 3 by those authors and our solution 1 was observed. Over a longer time interval, namely, 200 thousand years, the Earth’s rotation problem was solved initially towards the future [26]. The solution for, *ε* was found to feature a different oscillation structure and other onset times for extremes, and what was most important, the amplitude of new oscillations exceeding the amplitude of the previous oscillations by a factor of 7 - 8 was revealed. From this time on, a solution check for the Earth’s rotation problem was initiated. This check lasted for three years.

The Earth’s rotation problem is one of the most complex problems in mechanics. Its solution can depend on the fundamental assumptions made while deriving the equations, on the choice of initial conditions, and on the procedure of solution reduction to the Earth’s moving orbit. That is why a cardinal check of obtained results would be their obtaining, without solving differential equations for rotational motion.

While studying the orbits of the bodies, we found that the evolution of the Moon’s orbital axis was similar to that of Earth’s axis of rotation. That result had led us to a compound model for Earth rotation in which, part of the Earth’s mass was uniformly distributed among peripheral bodies rotating around a central body, along a circular orbit. Under the action of the Moon, the Sun, and the planets, the orbits of the peripheral bodies started changing. It should be noted that under the axis of the orbit is meant a perpendicular to its plane. The evolution of the orbital axis of one of the bodies, modeled the evolution of the Earth’s axis. Such modeling of the Earth’s rotational motion included performing several solution stages of the orbital problem being treated with the Galactica program. In the initial series of our studies [20 , 27], three models were studied, and the possibility of modeling the evolution of the Earth’s axis was confirmed. In those models, the precession periods of the orbit axes were 170 and 2604 years, whereas the average period of the Earth’s precession axis is known to be 25,740 years. Was subsequently yet developed 11 models, while has not been reached the required period of precession [7 , 13]. In the 13^{th} model, the orbital radius of the peripheral bodies was equal to the Earth’s radius, the rotation period of the bodies was 0.142 hour, and the interaction between the model’s bodies was amplified by a factor of 9.6 in comparison with the gravitational interaction. Thus, the bodies of the 13^{th} model rotated 170 times faster than the Earth does. For studying the evolution of such models, the Galactica program was developed; in this model, the possibility of changing the interaction between the model’s bodies was provided.

The solution of the problem for the 13^{th} compound model of the Earth, over a time interval of 300 years [7 , 13], has yielded all characteristics of the dynamics of the Earth’s axis, including the half-monthly and half-year oscillations of the angles *ε* and *ψ *and the oscillations of those angles with a 18.6-year period. The amplitudes of those oscillations have also proved to be coincident with the solution results of the Earth’s rotation problem. Such a coincidence in the results of the model problem with the results of the direct problem occurs once over a period of 3 thousand years. Further solution of the problem is hampered by the necessity to reduce the integration step to values for which the computing time turns out to be too long. Thus, over a time interval of 3000 years the compound model of the Earth has confirmed the obtained results of integration of the differential equation for Earth’s rotation. The latter indicates that the assumptions and simplifications adopted in the derivation of the equations, the derivation of the equations itself, the solution method, and the transforming of the integrated data into the final form have also been confirmed.

A second check consisted in using an alternative integration method. In the program DfEqAl1.for solving the Earth’s rotation problem, a fourth-order Runge-Kutta integration method as implemented by Krut’ko *et al.* [28] was used. Over a time interval of 200 thousand years, a growth of daily oscillations of the derivatives
$\stackrel{\dot{}}{\psi}$ and
$\stackrel{\dot{}}{\theta}$ was identified. Next, a code for solving the DfEqADP8.for program was developed, which uses the Dormand-Prince method, *i.e.* the eight-order Runge-Kutta integration method [29]. In integrating the equations of rotational motion over a period of 200 thousand years all previously obtained results have found confirmation. In the latter case, the amplitude of the daily oscillations of
$\stackrel{\dot{}}{\psi}$ and
$\stackrel{\dot{}}{\theta}$ showed no increase and remained at one and the same level. Thus, the particular method for integrating the equations does not affect the obtained results, and the application of a more accurate method confirms them.

A third check consisted in employing another method for solving the problem. The differential equations for rotational motion contained the coordinates of ten bodies acting on the Earth (eight planets, the Sun, and the Moon). Those coordinates were determined while solving the orbital problem with the Galactica program. However, the data array over large time intervals obtained with an integration step equal to 10^{−5} - 10^{−4} years took an overly large memory space. That is why we have developed a mathematical model for the Solar system [30] that yielded the coordinates of the bodies at a desired time, calculated by the formulas for elliptic motion in which the ellipse parameters at each time, were determined from the data initially calculated by the Galactica program. In the course of the solution process of the problem, the mathematical model for the Solar system was subjected to a thorough check. Nonetheless, there still existed a probability that, over large time intervals, the insignificant differences between the results of the mathematical model for the Solar system and the coordinate values obtained by the Galactica problem, could affect the evolution of the rotational parameters *ε* and *ψ*. We have developed a new program, glc3rte2.for, for joint solution of the orbital problem and the Earth’s rotation problem. In the latter program, in one step over time, the Galactica program solved the orbital problem and, then, the Dormand-Prince method was used to solve, in the same step, the Earth’s rotation problem. With the help of the new program, we have obtained solutions of those problems over different time intervals, including the interval of 200 thousand years. All previously obtained results have found confirmation. The latter check also confirms the mathematical model for the Solar system over large time intervals.

Table 1 gives a quantitative comparison of precession periods *P _{prN}*, and the minimum and maximum obliquities (

Table 1. Comparison of results obtained by three methods for integration of rotational-motion equations over a period of 200 thousand years: RK-4—Runge-Kutta method of the fourth order; DP-8— Runge-Kutta method of the eighth order in Dormand-Prince realization; Gal—the bodies’ coordinates are determined by the Galactica program, and the rotational-motion equations are solved by the DP-8 method.

third method, differed in terms of precession period *P _{prN}* in the fourth digit, and in terms of

Thus, the various tests and verification of the initial solution method of the Earth’s rotation problem and also, independent solutions of the same problem by three other methods, have confirmed the fact that the Earth’s rotational axis executes oscillations at an amplitude of 7 - 8 times greater than that obtained in the previous solutions.

4.4. Adequacy of the Solution of the Third and Forth Problems

As it was noted previously, the third problem on the determination of the obliquity *ε* between the moving equatorial plane 4 and the moving orbital plane 5, and the perihelion angle *φ _{pγ}* (Figure 2), was a geometrically complex problem because a multitude of revolutions on the axes
$\stackrel{\to}{N}$ and
$\stackrel{\to}{S}$

Those transformations were derived; yet, they involved inverse trigonometric functions, known to be many-valued. The expressions themselves are cumbersome, with some their parts presenting imaginary or infinitely growing expressions. Those singularities needed to be identified for the revealing factors underlying their behavior and for developing algorithms for elimination of these factors. Initially, this problem was solved by means of spherical geometry. However, because of the complexity of involved logical concepts, there was no firm belief in adequacy of the obtained solution. Fortunately, an idea suggesting a second method had emerged. The axial vectors, $\stackrel{\to}{N}$ and $\stackrel{\to}{S}$ were replaced with their projections onto the axes of a Cartesian coordinate system. Then, trigonometric means were employed to derive the required transformations. The two transformation systems have allowed us to reveal errors in each of the systems and fix them, unless both systems started yielding the same results in the examined 20-million-year time interval.

The new algorithm developed for calculating the insolation in the fourth problem was checked [19] by performing insolation calculations using the orbital data by Milankovitch [1], J. Laskar *et al.* [6]. The new algorithm has yielded results, the same as those yielded by the Milankovitch algorithm.

4.5. Physical Cause of the Difference between the New and Previous Theory

As it was noted previously, according to the new solutions for the angles *θ* and *ψ*, the Earth’s rotation axis
$\stackrel{\to}{N}$ rotates around the vector
${\stackrel{\to}{M}}_{2}$ (Figure 2), and it additionally exerts oscillations about that vector. We have calculated the angle *θ _{M}*

${P}_{rl}=\frac{{P}_{prS}\cdot {P}_{prN}}{{P}_{prS}-{P}_{prN}}$ . After substitution of actual values, we have obtained a value *P _{rl}* = −41.1 kyr, that is, the

very oscillation period of the angles *θ _{M}*

In the new theory (Figure 2), the vectors
${\stackrel{\to}{M}}_{2}$ and
$\stackrel{\to}{M}$ of the precessional axes have different orientations. That is why the moments of forces by which all bodies act on the Earth exhibit a wide range of oscillations. As a result, the oscillation range of *θ* with respect to the motionless orbit, (3) also increases. In addition, the oscillation range of the angle, *ε* between the moving orbital axis
$\stackrel{\to}{S}$ and the moving axis
$\stackrel{\to}{N}$ increases as well. All in all, in the new theory, the oscillation amplitude of angle, *ε* turns out to be exceeding the same amplitude in the previous theory by a factor of 7 - 8.

Thus, the oscillation period of obliquity *ε* in the previous theory of Earth-axis dynamics, was the period *P _{rl}* of the relative precession of Earth’s rotational axis
$\stackrel{\to}{N}$ and its orbit
$\stackrel{\to}{S}$ . That is why it is a fact that the precessional axes
$\stackrel{\to}{M}$ and
${\stackrel{\to}{M}}_{2}$ of the orbital axis
$\stackrel{\to}{S}$ and the Earth axis
$\stackrel{\to}{N}$ were assumed coincident, was the physical cause behind obtaining erroneous results in the previous theory.

4.6. Final Verification

The final check of the Astronomical theory of climate change consisted in the comparison of its results with paleoclimate. While analyzing data gained by geologists, geographers, and other specialists in the field of paleoclimate, we have found [22 , 31 , 32] four extremes of summer insolation ${Q}_{s}^{\text{65}N}$ having occurred over the past 50 thousand years; namely, 4.16, 15.88, 31.28, and 46.44 thousand years ago; those dates correspond to the middle of Holocene, to the middle of the last ice age, to the middle of a warm period, and to the middle of the penultimate ice maximum, respectively. Those events are called differently in different regions of the world; yet, all of them have left their traces in Siberia, Europe, and Northern America.

The whole complex of performed studies and their tests outlined here give us grounds to assert that, in the present article, results of Astronomical theory of climate change obtained by taking into account all studies performed during past centuries, are reported.

Changes of obliquity, *ε* and angle of the perihelion, *φ _{pγ}*, as well as the eccentricity of the orbit,

5. CONCLUSION AND FURTHER DEVELOPMENT OF THE ASTRONOMICAL THEORY OF CLIMATE CHANGE

As a result of the interaction of Solar-system bodies, evolution of the Earth’s orbital and rotational motions proceeds; this evolution in turn gives rise to insolation oscillations being the cause of climate changes observed over time intervals of tens of thousands of years. The same interactions lead to the evolution of the Sun’s motion about the center of mass of the Solar system [20 , 33] and also, to a change in the rotational motion of the Sun. Studies show [34], that the change in motions presents the cause of variations in solar activity. The radiation fluxes due to the Sun and its substance act on the Earth’s upper shells, thus leading to circulation changes arising in its atmosphere and ocean. Very likely, these factors act as the cause of short-term climate variations occurring over periods of several ten and hundred years. Further development of the Astronomical theory of climate change will be related to the determination of those oscillations.

ACKNOWLEDGEMENTS

The materials of this work were obtained as a result of research at the Institute of Earth’s Cryosphere, Tyum. SC of SB RAS, Federal Research Center for two decades, and in recent years, the research has been carried out under the theme 121041600047-2.

The results of the new Astronomical theory of climate change are based on the solution of the problems about the Earth’s orbital and rotational motion, obtained from supercomputers in the shared use center at the Siberian Supercomputer Center, Institute of Computational Mathematics and Mathematical Geophysics SB RAS, Novosibirsk, Russia.

English language within the paper has been edited by Walter Babin, researcher and founder of the General Science Journal, Canada.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

[1] | Milankovitch, M. (1930) Mathematische Klimalhre und Astronomische Theorie der Klimaschwankungen (Mathematical Climatology and the Astronomical Theory of Climate Change). Gebruder Borntraeger, Berlin. |

[2] | Brouwer, D. and Van Woerkom, A.J.J. (1950) The Secular Variation of the Orbital Elements of the Principal Planets. Astronomical Papers for American Ephemeris and Nautical Almanac, 13, 81-107. |

[3] | Sharaf, Sh.G. and Budnikova, N.A. (1969) The Secular Variation of the Orbital Elements of the Earth and the Astronomic Theory of Climate Change. Transactions, No. 14, 48-109. (In Russian) |

[4] |
Berger, A. and Loutre, M.F. (1991) Insolation Values for the Climate of the Last 10 Million Years. Quaternary Science Reviews, 1, 297-317. https://doi.org/10.1016/0277-3791(91)90033-Q |

[5] |
Edvardsson, S., Karlsson, K.G. and Engholm, M. (2002) Accurate Spin Axes and Solar System Dynamics: Climatic Variations for the Earth and Mars. Astronomy & Astrophysics, 384, 689-701. https://doi.org/10.1051/0004-6361:20020029 |

[6] |
Laskar, J., Robutel, P., Joutel, F., Gastineau, M., Correia, A.C.M. and Levrard, B. (2004) A Long-Term Numerical Solution for the Insolation Quantities of the Earth. Astronomy & Astrophysics, 428, 261-285. https://doi.org/10.1051/0004-6361:20041335 |

[7] |
Smulsky, J.J. (2014) Basic Positions and New Results of the Astronomical Theory of Climate Change. VINITI, Tyumen, No. 258-V2014, 30 p. (In Russian) http://www.ikz.ru/~smulski/Papers/OsPoATLP3.pdf |

[8] |
Smulsky, J.J. (2004) The Theory of Interaction. Cultural Information Bank, Ekaterinburg, 304 p. http://www.ikz.ru/~smulski/TVEnA5_2.pdf |

[9] |
Smulsky, J.J. (2012) Galactica Software for Solving Gravitational Interaction Problems. Applied Physics Research, 4, 110-123. https://doi.org/10.5539/apr.v4n2p110 |

[10] |
Smulsky, J.J. (2012) The System of Free Access Galactica to Compute Interactions of N-Bodies. International Journal of Modern Education and Computer Science, 11, 1-20. https://doi.org/10.5815/ijmecs.2012.11.01 |

[11] |
Smulsky, J.J. (2018) Future Space Problems and Their Solutions. Nova Science Publishers, New York, 269 p. http://www.ikz.ru/~smulski/Papers/InfFSPS.pdf |

[12] |
Smulsky, J.J. (2011) The Influence of the Planets, Sun and Moon on the Evolution of the Earth’s Axis. International Journal of Astronomy and Astrophysics, 1, 117-134. https://doi.org/10.4236/ijaa.2011.13017 |

[13] |
Smulsky, J.J. (2016) Fundamental Principles and Results of a New Astronomic Theory of Climate Change. Advances in Astrophysics, 1, 1-21. https://doi.org/10.22606/adap.2016.11001 |

[14] |
Smulsky, J.J. (2020) The Evolution of the Earth’s Rotational Movement for Millions Years. The Complex Systems, 1, 3-42. http://www.ikz.ru/~smulski/Papers/EVDZ03_1EnJc.pdf |

[15] |
Smulsky, J.J. and Krotov, O.I. (2014) New Computing Algorithm of the Earth’s Insolation. Applied Physics Research, 6, 56-82. https://doi.org/10.5539/apr.v6n4p56 |

[16] |
Newcomb, S. (1895) The Elements of the Four Inner Planets and the Fundamental Constants of Astronomy. Government Printing Office, Washington DC, 202 p. https://doi.org/10.1086/102176 |

[17] | Simon, J.L., Bretagnon, P., Chapront, J., Chapront-Touze, M., Francou, G. and Laskar, J. (1994) Numerical Expression for Precession Formulae and Mean Elements for the Moon and the Planets. Astronomy & Astrophysics, 282, 663-683. |

[18] |
Smulsky, J.J. (2018) The Sun’s Movement in the Sky Now and in the Past. Open Access Library Journal, 5, e4250. https://doi.org/10.4236/oalib.1104250 |

[19] |
Smulsky, J.J. (2018) New Astronomical Theory of Ice Ages. LAP Lambert Academic Publishing, Riga, 132 p. (In Russian) http://www.ikz.ru/~smulski/Papers/InfNwATLPEn.pdf |

[20] |
Melnikov, V.P. and Smulsky, J.J. (2009) Astronomical Theory of Ice Ages: New Approximations. Solutions and Challenges. Academic Publishing House, Novosibirsk. http://www.ikz.ru/~smulski/Papers/AsThAnE.pdf |

[21] |
Smulsky, J.J. (2016) Evolution of the Earth’s Axis and Paleoclimate for 200,000 Years. LAP Lambert Academic Publishing, Saarbrucken, 228 p. (In Russian) http://www.ikz.ru/~smulski/Papers/InfEvEAPC02MEn.pdf |

[22] |
Smulsky, J.J. (2016) New Results on the Earth Insolation and Their Correlation with the Late Pleistocene Paleoclimate of West Siberia. Russian Geology and Geophysics, 57, 1099-1110. https://doi.org/10.1016/j.rgg.2016.06.009 |

[23] |
Grebenikov, E.A. and Smulsky, J.J. (2007) Evolution of the Mars Orbit on Time Span in Hundred Millions Years. Reports on Applied Mathematics. Russian Academy of Sciences: A. A. Dorodnitsyn Computing Center, Moscow, 63 p. (In Russian) http://www.ikz.ru/~smulski/Papers/EvMa100m4t2.pdf |

[24] |
Laskar, J. (1996) Marginal Stability and Chaos in the Solar System. In: Ferraz Mello, S., et al., Eds., Dynamics, Ephemerides and Astrometry of the Solar System, IAU, Amsterdam, 75-88. https://doi.org/10.1017/S0074180900127160 |

[25] | Bretagnon, P., Rocher, P. and Simon, J.L. (1997) Theory of the Rotation of the Rigid Earth. Astronomy and Astrophysics, 319, 305-317. |

[26] |
Smul’skii, I.I. (2013) Analyzing the Lessons of the Development of the Orbital Theory of the Paleoclimate. Herald of the Russian Academy of Sciences, 83, 46-54. https://doi.org/10.1134/S1019331613010073 |

[27] |
Mel’nikov, V.P., Smul’skii, I.I. and Smul’skii, Ya.I. (2008) Compound Modeling of Earth Rotation and Possible Implications for Interaction of Continents. Russian Geology and Geophysics, 49, 851-858. https://doi.org/10.1016/j.rgg.2008.04.003 |

[28] | Krutko, P.D., Maksimov, A.I. and Skvortsov, L.M. (1988) Automated System Engineering Algorithms and Programs. Radio and Communication, Moscow, 298 p. (In Russian). |

[29] |
Hairer, E., Nersett, S. and Vanner, G. (1990) Solving Ordinary Differential Equations. Mir, Moscow, 513 p. (In Russian) https://doi.org/10.1007/978-3-662-09947-6 |

[30] |
Smulsky, J.J. (2007) A Mathematical Model of the Solar System. In: Theoretical and Applied Problems in Nonlinear Analysis, Russian Academy of Sciences: A.A. Dorodnitsyn Computing Center, Moscow, 119-138. (In Russian) http://www.ikz.ru/~smulski/Papers/MatMdSS5.pdf |

[31] |
Smulsky, J.J. and Ivanova, A.A. (2018) Experience of Paleoclimate Reconstruction on Insolation Change on Example of Western Siberia in the Late Pleistocene. Climate and Nature, 1, 3-21. (In Russian) http://www.ikz.ru/~smulski/Papers/OpRcnPClmt6J.pdf |

[32] |
Smulsky, J.J. (2019) The Upcoming Tasks of Fundamental Science. Sputnik + Publishing House, Moscow, 134 p. (In Russian) http://www.ikz.ru/~smulski/Papers/InfPrZaFN.pdf |

[33] |
Smulsky, J.J. (2017) Cosmic Impacts on the Earth and Their Influence on the Arctic. The Complex Systems, 4, 27-42. (In Russian) http://www.ikz.ru/~smulski/Papers/CsmAcEIA.pdf |

[34] | Mörner, N.-A. (2016) Planetary Influence on the Sun and the Earth, and a Modern Book-Burning. Nova Science Publishers, New York, 196 p. |

Journals Menu

Contact us

customer@scirp.org | |

+86 18163351462(WhatsApp) | |

1655362766 | |

Paper Publishing WeChat |

Copyright © 2023 by authors and Scientific Research Publishing Inc.

This work and the related PDF file are licensed under a Creative Commons Attribution 4.0 International License.