© ESO, 2012

1. Introduction

Pierre Simon Laplace (1749–1827) in his Oeuvres Complètes (De Laplace 1878) was one of the first scientists to mention the notion of an invariable plane of the solar system. He introduced what seemed to be a natural reference plane when studying the motion of celestial bodies (comets, etc.). He showed that for any isolated N-body system, subject to mutual gravitational interactions, one can define an invariable plane that obeys a simple geometrical property and is the consequence of the dynamics of a conservative system. Thus, it is the natural reference plane when studying the system.

Over the past two centuries, a few astronomers (see for instance, See 1904; Innes 1920; Clemence & Brouwer 1955) have worked on determining the orientation of the invariable plane. In (See 1904), prior to the discovery of Pluto, the position of this plane was given with respect to the ecliptic of the epoch 1850.0. Its inclination was set to ${\mathrm{1}}^{\mathrm{◦}}{\mathrm{35}}^{\mathrm{\prime }}\mathrm{7}\stackrel{\mathrm{″}}{.}\mathrm{745}$ and the longitude of its ascending node to ${\mathrm{106}}^{\mathrm{◦}}{\mathrm{8}}^{\mathrm{\prime }}\mathrm{46}\stackrel{\mathrm{″}}{.}\mathrm{688}$. In Innes (1920), also prior to the discovery of Pluto, the orientation of the invariable plane was given with respect to both the ecliptic and the equator of the epoch 1900.0. The inclination and the longitude of the ascending node of the invariable plane with respect to the ecliptic (the equator, respectively) are ${\mathrm{1}}^{\mathrm{◦}}{\mathrm{34}}^{\mathrm{\prime }}\mathrm{59}\stackrel{\mathrm{″}}{.}\mathrm{42}$ and ${\mathrm{106}}^{\mathrm{◦}}{\mathrm{35}}^{\mathrm{\prime }}\mathrm{1}\stackrel{\mathrm{″}}{.}\mathrm{08}$ (${\mathrm{23}}^{\mathrm{◦}}{\mathrm{2}}^{\mathrm{\prime }}\mathrm{51}\stackrel{\mathrm{″}}{.}\mathrm{40}$ and ${\mathrm{3}}^{\mathrm{◦}}{\mathrm{52}}^{\mathrm{\prime }}\mathrm{41}\stackrel{\mathrm{″}}{.}\mathrm{38}$, respectively).

The last known paper investigating the topic is (Burkhardt 1982). In this paper, the orientation of the invariable plane was computed using the numerical Development Ephemerides DE96, DE102, DE108, and DE110. The inclination and the longitude of the ascending node of the invariable plane were determined with respect to both the ecliptic-equinox and the equator-equinox of the epoch B1950.0. The orientation was also given with respect to the equator-equinox and ecliptic-equinox of the epoch J2000.0, by applying the precession matrix given in Lieske et al. (1977) and Lieske (1979). For other papers investigating the characteristics of the invariable plane with respect to the mean ecliptic-equinox of the epoch B1950.0, we refer to (Burkhardt 1982). These papers underline the principal difficulty in determining its orientation, which lies in the uncertainties of the masses of all the objects involved as well as in their positions and velocities.

The aim of this paper is to deepen the investigation of the notion of the invariable plane of the solar system and to define as accurately as possible its orientation with respect to the ecliptic and the equator of the epoch J2000.0. There have been no recent publications investigating the concept since (Burkhardt 1982), though it has been used recently, in the construction of the INPOP ephemeris. Indeed, the ring modeling the perturbation of the main-belt asteroids is assumed to be located in the invariable plane of the ephemeris (Kuchynka et al. 2010; Kuchynka 2010).

In this paper, we propose the use of the most recent long-term numerical ephemeris data to determine the orientation of the invariable plane. In addition to the Sun and its eight planets (the Earth being replaced by the Earth-Moon barycentre – EMB), we show that it is necessary to take account of the motions of dwarf planets such as (134340) Pluto and (1) Ceres, as well as the asteroids (2) Pallas and (4) Vesta in order to optimise this determination. In this paper, we also estimate the individual contributions of each body constituting the system. Our study provides a non-negligible update on the topic as well as a non-standard test of the numerical ephemerides. We note that the invariable plane is examined by the IAU/IAG Working Group on Cartographic Coordinates and Rotational Elements in 2006 (Seidelmann et al. 2007) and in the report of WGCCRE 2009 (Archinal et al. 2011) as a reference plane to determine the obliquity of solar system objects. For the planets, the primary pole is the one on the north side of the invariable plane of the solar system.

2. Invariable plane

2.1. Definition of the invariable plane

Considering the solar system as isolated, its total angular momentum vector is constant with respect to both spatial and time coordinates. Thus, the invariable plane¹ is defined as the plane perpendicular to the total angular momentum vector of the solar system that passes through its barycentre. Being fixed, it provides a permanent natural reference plane, whereas the ecliptic slightly moves with time.

In this paper, we determine the orientation of the invariable plane by setting its inclination and the longitude of its ascending node with respect to both the ICRF (origin and equator) and the ecliptic-equinox of the epoch J2000.0.

2.2. Determining the orientation of the invariable plane

In Newtonian mechanics, the total angular momentum vector of an N-body system (here, we disregard the rotation of all the bodies) is given by $\begin{array}{ccc}{L}\mathrm{tot}\mathrm{=}\sum _{\mathit{j}\mathrm{=}\mathrm{1}}^{\mathit{N}}{\mathit{m}}_{\mathit{j}}{r}\mathit{j}\mathrm{\times }{ṙ}\mathit{j}\mathit{,}& & \end{array}$(1)where m_j, r_j, and ṙ_j are the mass, barycentric position vector, and barycentric velocity vector of the jth body, respectively. The norm of the total angular momentum vector is given by $\begin{array}{ccc}{\mathit{L}}_{\mathrm{tot}}\mathrm{=}\sqrt{{{\mathit{L}}^{\mathrm{2}}}_{\mathrm{1}}\mathrm{+}{{\mathit{L}}^{\mathrm{2}}}_{\mathrm{2}}\mathrm{+}{{\mathit{L}}^{\mathrm{2}}}_{\mathrm{3}}}\mathit{,}& & \end{array}$(2)with $\begin{array}{c}\\ {\mathit{L}}_{\mathrm{1}}\mathrm{=}\sum _{\mathit{j}\mathrm{=}\mathrm{1}}^{\mathit{N}}{\mathit{m}}_{\mathit{j}}\mathrm{(}{\mathit{y}}_{\mathit{j}}\mathit{ż}\mathit{j}\mathrm{-}{\mathit{z}}_{\mathit{j}}\mathit{ẏ}\mathit{j}\mathrm{)}\mathit{,}\\ {\mathit{L}}_{\mathrm{2}}\mathrm{=}\sum _{\mathit{j}\mathrm{=}\mathrm{1}}^{\mathit{N}}{\mathit{m}}_{\mathit{j}}\mathrm{(}{\mathit{z}}_{\mathit{j}}\mathit{ẋ}\mathit{j}\mathrm{-}{\mathit{x}}_{\mathit{j}}\mathit{ż}\mathit{j}\mathrm{)}\mathit{,}\\ {\mathit{L}}_{\mathrm{3}}\mathrm{=}\sum _{\mathit{j}\mathrm{=}\mathrm{1}}^{\mathit{N}}{\mathit{m}}_{\mathit{j}}\mathrm{(}{\mathit{x}}_{\mathit{j}}\mathit{ẏ}\mathit{j}\mathrm{-}{\mathit{y}}_{\mathit{j}}\mathit{ẋ}\mathit{j}\mathrm{)}\mathit{,}\end{array}$(3)where x_j, y_j, z_j (ẋ_j, ẏ_j, ż_j, respectively) are the components of the barycentric position (velocity, respectively) vector of the jth body.

From Brouwer & Clemence (1961, Chap. I), we have $\begin{array}{c}\\ {\mathit{L}}_{\mathrm{1}}\mathrm{=}{L}\mathrm{tot}\hspace{0.17em}\sin \mathrm{\Omega }\hspace{0.17em}\sin \mathit{i,}\\ {\mathit{L}}_{\mathrm{2}}\mathrm{=}\mathrm{-}{L}\mathrm{tot}\hspace{0.17em}\cos \mathrm{\Omega }\hspace{0.17em}\sin \mathit{i,}\\ {\mathit{L}}_{\mathrm{3}}\mathrm{=}{L}\mathrm{tot}\hspace{0.17em}\cos \mathit{i,}\end{array}$(4)where L₁, L₂, L₃ are the cartesian barycentric components of the total barycentric angular momentum vector L_tot given by (1) and Ω, i are, respectively, the longitude of the ascending node and the inclination of the invariable plane with respect to the reference frame in which are given the position and velocity vectors of the N bodies of the system. In our case, the positions and velocities vectors provided by the numerical ephemerides, are given with respect to the ICRF.

Using Eqs. (1) to (4), we compute the values of i$\begin{array}{c}\\ \mathit{i}\mathrm{=}\arccos \frac{{\mathit{L}}_{\mathrm{3}}}{{\mathit{L}}_{\mathrm{tot}}}\end{array}$(5)and Ω $\begin{array}{c}\\ \mathrm{\Omega }\mathrm{=}\arctan \frac{\mathrm{-}{\mathit{L}}_{\mathrm{1}}}{{\mathit{L}}_{\mathrm{2}}}\mathrm{·}\end{array}$(6)The relativistic effects being taken into account up to first order in the post-Newtonian approximation in the compilation of the numerical ephemerides DE405/DE406 (Standish 1998) and INPOP10a (Fienga et al. 2010), (Kuchynka 2010), the total angular momentum vector initially given by (1) becomes $\begin{array}{ccc}{L}\mathrm{tot}\mathrm{=}\sum _{\mathit{j}\mathrm{=}\mathrm{1}}^{\mathit{N}}{\mathit{m}}_{\mathit{j}}^{\mathrm{\ast }}{r}\mathit{j}\mathrm{\times }{ṙ}\mathit{j}\mathit{,}& & \end{array}$(7)where m_j is the mass of the jth body being replaced by a corresponding effective mass ${\mathit{m}}_{\mathit{j}}^{\mathrm{\ast }}$ (the so-called Tolmann-Mass) $\begin{array}{ccc}{\mathit{m}}_{\mathit{j}}^{\mathrm{\ast }}\mathrm{=}{\mathit{m}}_{\mathit{j}}\mathit{.}\left[\mathrm{1}\mathrm{+}\frac{\dot{{r}\begin{array}{}\mathrm{2}\\ \mathit{j}\end{array}}}{\mathrm{2}{\mathit{c}}^{\mathrm{2}}}\mathrm{-}\frac{\mathrm{1}}{\mathrm{2}{\mathit{c}}^{\mathrm{2}}}\left(\sum _{\mathit{k}\ne \mathit{j}}\frac{\mathit{G}{\mathit{m}}_{\mathit{k}}}{\mathrm{|}{r}\mathit{k}\mathrm{-}{r}\mathit{j}\mathrm{|}}\right)\right]\mathit{,}& & \end{array}$(8)c = 2.99792458 × 10⁸ m s^-1 is the velocity of light in vacuum (Luzum et al. 2011) and G is the gravitational constant (see Standish et al. 1976, for more details about (7); a proof can also be found in Brumberg 1991).

2.3. Method

To determine the position of the invariable plane, we use two different long-term numerical ephemerides:

• The Development Ephemerides DE 405/406(Standish 1998).

• The INPOP10a Ephemeris (Fienga et al. 2010; Kuchynka 2010).

The Development Ephemerides DE405/DE406 (Standish 1998) are not the most recent ephemerides provided by the JPL, but are the ones with the longer available time span (about 6000 y for DE406). The orientation of the inner planetary system of DE405 is accurate to about 0.001 arcsec. The position ephemerides of the outer planets provided by the DE405/406 rely mainly upon optical observations. For instance, Jupiter’s observations are fit by DE405 to an accuracy of 0.01 arcsec. The perturbations of the 300 asteroids that have the strongest influence on the Earth, Moon, and Mars were introduced, in two steps:

• 1.

  The three largest main belt asteroid objects (1) Ceres,(2) Pallas, and (4) Vesta were integrated as massiveparticles under the gravitational perturbation of the Sun, the eightplanets, the Moon, and Pluto. The results were fit by temporaryChebyshev polynomials.

• 2.

  The 297 remaining main belt asteroids, considered as massless, were integrated under the gravitational influence of the Sun, the eight planets, the Moon, Pluto, Ceres, Pallas, and Vesta. The final results were then fit with Chebyshev polynomials, providing us with a full planetary and lunar integration.

INPOP10a (Fienga et al. 2010; Kuchynka 2010) is the most recent available numerical ephemeris provided by the IMCCE. Thus, it exploits the most recent observational constraints provided by the MEX and VEX missions for Mars and Venus, the Cassini mission for Saturn, several sets of flybys data of Jupiter, Uranus, and Neptune (Pionner 10 and 11, Viking 1 and 2, Ulysses and Cassini), and stellar occultations in the case of Pluto. INPOP10a is fitted over the time interval  [1914.2:2009.7] , the results are then extrapolated over a longer time interval (Table 1). In addition to the 300 most strongly perturbing main belt asteroids, the perturbations of 24635 asteroids are averaged over a ring, of an analytically estimated radius of 2.8 AU (Kuchynka et al. 2010; Kuchynka 2010). The initial inclination of this ring, with respect to the invariable plane³, was chosen to be equal to zero.

Table 1 illustrates the available time span for each numerical ephemeris. Moreover, the primary constants for both the DE(405/406) and INPOP10a ephemerides are given in Table 2.

Using Eqs. (2) to (8), we computed from each ephemeris over its entire time span, with a 1 d time step, the values of i and Ω (the inclination (Eq. (5)) and the longitude of the ascending node (Eq. (6)), respectively) of the invariable plane with respect to the origin-equator of the ICRF. For this purpose, we considered two physical systems, one basic and the other more complete. This is described in detail in Sect. 3.

The ephemerides being given with respect to the ICRF, the resulting Ω and i of the invariable plane give its orientation with respect to this reference frame, where the origin of the longitude of the ascending node is the ICRF’s one (intersection of the orbit with the x-axis of the ICRF, Fig. 1), and the inclination is given with respect to the equator of the ICRF given by the x − y axes.

To obtain the position of the invariable plane with respect to the equinox-ecliptic of the epoch J2000.0, we proceed to a coordinate transformation for the total angular momentum vector L (see Fig. 1).

Noting L_ecl (L_ICRF, respectively), the total angular momentum of the solar system computed with respect to the equinox-ecliptic J2000.0 (the ICRF, respectively), this transformation involving two matrices is given by the equation $\begin{array}{ccc}{L}\mathrm{ecl}\mathrm{=}{\mathit{R}}_{\mathit{x}}\mathrm{\left(}\mathit{ϵ}\mathrm{\right)}{\mathit{R}}_{\mathit{z}}\mathrm{\left(}\mathit{\varphi }\mathrm{\right)}{L}\mathrm{ICRF}& & \end{array}$(9)with: $\begin{array}{ccc}{\mathit{R}}_{\mathit{z}}\mathrm{\left(}\mathit{\varphi }\mathrm{\right)}\mathrm{=}\left(\begin{array}{c}\\ & & \\ & & \\ & & \end{array}\right)& & \end{array}$and $\begin{array}{ccc}{\mathit{R}}_{\mathit{x}}\mathrm{\left(}\mathit{ϵ}\mathrm{\right)}\mathrm{=}\left(\begin{array}{c}\\ & & \\ & & \\ & & \end{array}\right)\mathit{,}& & \end{array}$where ϵ is the obliquity, that is to say the angle between the Earth’s equator and the ecliptic also known as the Earth’s axial tilt, and ϕ is a longitude bias, since the origin of the ICRF does not correspond to the equinox of the epoch J2000.0.

The values of the parameters ϵ and ϕ, for the DE405/406 and INPOP10a ephemerides are given in Table 3 (Simon 2011, priv. comm.).

3. Results

For each ephemeris data set, we computed the inclination i and the longitude of the ascending node Ω of the invariable plane, with a time step of 1d, with respect to both:

• the ICRF, by using Eqs. (1) to (8);

• the equinox-ecliptic of the epoch J2000.0, after applying the rotation given by Eq. (9).

We have distinguished two cases:

• A basic physical system including the Sun, the dwarf planet(134340) Pluto, and the eight planets (Souami &Souchay 2011). This system was the oneconsidered by (Burkhardt 1982).

• The previous system, to which we added the dwarf planet (1) Ceres as well as the two asteroids (4) Vesta and (2) Pallas. This allows us to evaluate the effects of these three additional small bodies.

3.1. The basic system: the Sun, the eight planets, and Pluto

Our N-body system (N = 10) is defined as follows: the Sun, (134340) Pluto, and the eight planets, the Earth being replaced by the EMB. This N-body system is identical to the one considered in previous studies, hence, we investigate it and compare our results to (Burkhardt 1982).

3.1.1. Results at epoch J2000.0

In Table 4, we present our two sets of measurements of the orientation of the invariable plane at the epoch J2000.0, obtained by using both the DE405/DE406 and INPOP10a ephemerides, finding good agreement between thems. In particular, when comparing the results obtained by using DE405 and INPOP10a, we find an agreement of up to 3 mas (milliarcseconds) for i and 37 mas for Ω with respect to the ICRF (origin-equator). The differences are probably due to the different values of the primary parameters used in the two numerical ephemerides (see Table 2).

As for the differences found between the results observed with respect to the ecliptic-equinox of J2000.0, the differences that we found here were amplified by the change in the determination of the parameters ϕ and ϵ. They are respectively 12 mas for i and 295 mas for Ω.

When comparing the results obtained by using the DE405 and DE406 ephemerides, we observe an agreement between them up to 0.1 mas and 0.01 mas for i with respect to both the ICRF and the ecliptic-equinox of the epoch J2000.0. For Ω, the agreement is as close as 0.02 mas and 0.16 mas.

The last two columns in Table 4 show the results obtained by (Burkhardt 1982), given with respect to the equator-equinox system of J2000.0, by using the precession matrix between B1950.0 and J2000.0. When comparing these results to ours, obtained with respect to the ICRF, we get large differences of $\mathrm{9}\stackrel{\mathrm{″}}{.}\mathrm{8}$ for i, and 1′14′′ for Ω. This can still be explained by the use of different ephemerides data sets with different values of the primary constants, but the large part of the difference is undoubtedly due to the use of a different reference plane used in the computations. In other words, the origin and equator of the ICRF do not coincide exactly with the equinox-equator.

3.1.2. Long-term accuracy of the determination of the invariable plane

The position and velocity vectors of each body, provided by the adopted numerical ephemerides, are subject to a limited precision. Moreover, neglecting bodies such as the largest asteroids induces some variations in the orientation of the invariable plane with respect to time, when ideally it should remain fixed as stated in Sect. 2. By studying the temporal variations in the invariable plane’s orientation, we check the accuracy of its determination in the case of the basic system, with respect to both the ICRF and the ecliptic-equinox of the epoch J2000.0, over the entire available time span (Table 1) of each of the numerical ephemerides DE405, DE406, and INPOP10a, using a time step of 1d. We show the results in Table 5. In particular, we give the maximum and minimum values of i and Ω, with respect to both the ICRF and the equinox-ecliptic of the epoch J2000.0. We also provide the value of the maximal variations in Δi and ΔΩ in all cases. Here again, our results show a good agreement between INPOP10a and DE405/406.

Moreover, the amplitude of the temporal variations Δi (respectively ΔΩ) does not exceed 1.22 mas (503 mas, respectively) with respect to the ecliptic-equinox J2000.0. The maximum amplitude of these variations with respect to the ICRF does not exceed 13.59 mas (7.53 mas, respectively).

With respect to the ICRF, we clearly observe a linear trend of i as a function of time, at a rate of $\mathrm{-}\mathrm{2}\stackrel{\mathrm{″}}{.}\mathrm{365252}\mathrm{\times }{\mathrm{10}}^{-6}\mathit{/}\mathit{y}$ ($\mathrm{-}\mathrm{2}\stackrel{\mathrm{″}}{.}\mathrm{255406}\mathrm{\times }{\mathrm{10}}^{-6}\mathit{/}\mathit{y}$, respectively) over the entire available time span of the ephemeris DE405 (Fig. 2(a)) and INPOP10a (Fig. 2(b)). Nevertheless, in the second case, a quadratic component is necessary to determine the closest fit to the curve. We note that the two slopes are of a negative sign and values very close to each other.

These linear behaviours can still be discerned with respect to the ecliptic-equinox of J2000.0 (Figs. 3(a), b), with slopes of positive sign and values close to each other, namerly $\mathrm{8}\stackrel{\mathrm{″}}{.}\mathrm{4037}\mathrm{\times }{\mathrm{10}}^{-5}\mathit{/}\mathit{y}$ for the DE405 ephemeris (Fig. 3(a)) and $\mathrm{8}\stackrel{\mathrm{″}}{.}\mathrm{9197}\mathrm{\times }{\mathrm{10}}^{-5}\mathit{/}\mathit{y}$ for the INPOP10a ephemeris (Fig. 3(b)).

The temporal variations in Ω determined from the DE405 ephemeris over its 600 y time span, with respect to the origin-equator of the ICRF and the ecliptic-equinox of J2000.0 are plotted in Figs. 2(c), 3(c), respectively. In the first case, we observe that, the linear slope has a much lower absolute value ($\mathrm{1}\stackrel{\mathrm{″}}{.}\mathrm{4096694}\mathrm{\times }{\mathrm{10}}^{-6}\mathit{/}\mathit{y}$ ) with an opposite sign, than in the second case ($\mathrm{8}\stackrel{\mathrm{″}}{.}\mathrm{40370}\mathrm{\times }{\mathrm{10}}^{-5}\mathit{/}\mathit{y}$).

Moreover, when studying these temporal variations by using INPOP10a, we find that the signals are clearly dominated by a quadratic drift with respect to both the origin-equator of the ICRF (Fig. 2(d)) and the ecliptic-equinox of J2000.0 (Fig. 3(d)).

Finally, we study Δi and ΔΩ, when using the long-term Development Ephemeris DE406 (Standish 1998) over its 6000-year time span (Table 1), with respect to the ICRF. Results are plotted in Figs. 4 and 5, respectively. A linear behaviour is observed in both cases. With respect to the ICRF, the rate of Δi is evaluated to be $\mathrm{-}\mathrm{2}\stackrel{\mathrm{″}}{.}\mathrm{265424036}\mathrm{\times }{\mathrm{10}}^{-6}$/yr, which is close to the value found above with the INPOP10a ephemeris over its 2053 y time interval.

We confirm the linear behaviour of ΔΩ, already observed in Fig. 2(c) for DE405, but do not indentify a quadratic component similar to that found for INPOP10a Fig. 2(d).

We note that (Burkhardt 1982), using the DE102 ephemeris over the time interval [1497, 2249], also identified a linear behaviour of the temporal variations of i and Ω of the invariable plane.

3.1.3. Short-term accuracy of the determination

In addition to the linear and quadratic components, we note in Figs. 2(a)–d and 3(a)–d the small temporal variations, at high frequencies, in the orientation of the invariable plane.

In Figs. 6 and 7, we show these variations over the 100 y time interval [1950.0, 2050.0], with respect to the ICRF. We note the very good agreement between the results obtained by using the DE405 and the INPOP10a ephemerides, for both Δi (Fig. 6) and ΔΩ (Fig. 7).

Since the variations in i and Ω have the same profile, even when derived from two independent ephemerides, we deduce that they are probably induced by an imperfection in the computation of the total angular momentum, rather than in the ephemerides themselves. This is investigated in Sect.3.2.

3.2. A more complete system

If our ten-body system were an isolated system, its total momentum would be conserved and the values of the longitude of the ascending node Ω and the inclination i of the invariable plane would be constant over time. This is not the case as can be seen in Figs. 2(a)–d with respect to the ICRF, and Figs. 3(a)–d with respect to the ecliptic-equinox of the epoch J2000.0. This is partially due to the uncertainties in the positions, velocities, and the masses of these bodies. This is also probably because, by considering only the planets, we have ignored the effect of smaller bodies, in particular the larger ones of the asteroid main belt. As previously stated in Sect. 2.3, (1) Ceres, (2) Pallas, and (4) Vesta were taken into account, being massive bodies, in the compilation of the numerical ephemerides DE405/406 (Standish 1998) and INPOP10a (Fienga et al. 2010; Kuchynka 2010). Thus, they should also be included in the computation of the total angular momentum of the system. The effect of these three bodies on the orientation of the invariable plane is investigated here.

For this purpose, we add to our previous basic system the three largest bodies of the asteroid main belt mentioned above and assess their effect on the orientation of the invariable plane. Thus, our 10-body system defined in Sect.3.1 becomes the following 13-body system: the Sun, the eight planets (the Earth being replaced by the Earth-Moon barycentre), the two dwarf planets (134340) Pluto and (1) Ceres, and the two asteroids (2) Pallas and (4) Vesta. Nevertheless, their positions and velocities are not directly available from the ephemerides.

Therefore, unlike the planetary ephemerides that were retrieved through the Access Direct files of Chebyshev polynomial fits to the cartesian positions and velocities of the planets, the Sun and the Moon, those of (1) Ceres, (2) Pallas, and (4) Vesta were collected from the IMCCE Ephemeris Web Server (IMCCE 2011).

This web server provides position ephemerides for the small bodies of the solar system by numerical integration, using the numerical planetary ephemerides mentioned above for the initial conditions of the massive bodies: the Sun, the planets (the Earth being replaced by the EMB), as well as the dwarf planet (134340) Pluto.

3.2.1. Results at the epoch J2000.0

Here, as in Sect. 3.1.1, we present the values of the orientation parameters i and Ω of the invariable plane at the epoch J2000.0, by taking into account the additional effects of (1) Ceres, (2) Pallas, and (4) Vesta, thus in the context of the more complete system.

The results obtained by using the DE405 and the INPOP10a ephemerides, with respect to both the origin-equator of the ICRF and the ecliptic-equinox of the epoch J2000.0 are given in Table 7.

Results with respect to the ICRF show the close agreement between the two ephemerides at the level of 4 mas for i, and 37 mas for Ω. As stated in Sect. 3.1.1, different primary body masses (see Table 2) were used in the compilation of the ephemerides. For instance, the masses of (1) Ceres, (2) Pallas, and (4) Vesta in 10^-10 × M_⊙ unit are given to the first decimal in the case of DE405 ephemeris, whereas they are given to the fifth decimal in the case of the INPOP10a. This might explain the small differences found in the orientation, between the two ephemerides.

With respect to the ecliptic-equinox of the epoch J2000.0, these differences are larger: 12 mas for i and 314 mas for Ω. This was already the case for the basic planetary system seen in Sect. 3.1.1. The ICRF (origin-equator) to equinox-ecliptic transformation in Eq. (9) with the parameters given in Table 3 (Simon 2011, priv. comm.), is most likely the origin of these amplifications. We note that the four values above showing the differences between the DE405 and INPOP10a results, are quasi-identical to those found with the basic system (Sect. 3.1.1).

Finally, we compare the results of the orientation of the invariable plane obtained at the epoch J2000.0, when considering the basic system (Table 4), to those obtained when considering the more complete system (Table 7). The effect of the three additional bodies is of the order of 4 mas on i and 22 mas on Ω with respect to the ICRF, for both the DE405 and INPOP10a ephemerides; with respect to the equinox-ecliptic of the epoch J2000.0, this effect is instead found to be 9 mas for i for both ephemerides, and 43 mas (62 mas) for Ω with DE405 (INPOP10a, respectively).

3.2.2. Results over 100 years

We have studied in Sect. 3.1.3 small temporal variations in the orientation of the invariable plane over the 100 y time interval [1950.0, 2050.0], in the case of the basic system.

The purpose of this subsection is to study the specific influence of the three largest Main Belt Asteroid objects, over the 100 y considered, and evaluate the possible improvements in the conservation of the angular momentum vector when taking them into account.

We introduce (1) Ceres, (4) Vesta, and (2) Pallas, one after the other in this order, the most massive to the least massive (see Table 2), into our previous 10-body system. In each case, we evaluate the peak-to-peak variations in both i and Ω of the invariable plane, with respect to the ICRF, using both the DE405 and INPOP10a ephemerides, over the time interval [1950.0, 2050.0]. We report these values in Table 10.

In the case of the basic system, these variations are $\mathrm{2}\stackrel{\mathrm{″}}{.}\mathrm{14}\mathrm{\times }{\mathrm{10}}^{-4}$ and $\mathrm{3}\stackrel{\mathrm{″}}{.}\mathrm{77}\mathrm{\times }{\mathrm{10}}^{-4}$ for i in the case of DE405 and INPOP10a, respectively. After introducing the dwarf planet (1) Ceres we observe, over the same 100 y time interval, a considerable reduction to $\mathrm{0}\stackrel{\mathrm{″}}{.}\mathrm{75}\mathrm{\times }{\mathrm{10}}^{-4}$ and $\mathrm{0}\stackrel{\mathrm{″}}{.}\mathrm{84}\mathrm{\times }{\mathrm{10}}^{-4}$, that is to say an improvement by a factor of 2.85 for DE405 and 4.48 for INPOP10a.

In contrast, for Ω, we note a small increase in the variations for both DE405 and INPOP10a. This is shown in Figs. 9(a) and b, respectively.

Nevertheless, after introducing the asteroid (4) Vesta in addition to (1) Ceres (Table 10), we note a significant improvement for Ω. In the case of the DE405 (INPOP10a, respectively) ephemeris, the peak-to-peak value goes from $\mathrm{2}\stackrel{\mathrm{″}}{.}\mathrm{33}\mathrm{\times }{\mathrm{10}}^{-4}$ ($\mathrm{2}\stackrel{\mathrm{″}}{.}\mathrm{43}\mathrm{\times }{\mathrm{10}}^{-4}$, respectively ) for the basic system, to $\mathrm{3}\stackrel{\mathrm{″}}{.}\mathrm{44}\mathrm{\times }{\mathrm{10}}^{-4}$ ($\mathrm{5}\stackrel{\mathrm{″}}{.}\mathrm{15}\mathrm{\times }{\mathrm{10}}^{-4}$, respectively) after introducing only (1) Ceres and decreases to $\mathrm{1}\stackrel{\mathrm{″}}{.}\mathrm{10}\mathrm{\times }{\mathrm{10}}^{-4}$ ($\mathrm{2}\stackrel{\mathrm{″}}{.}\mathrm{41}\mathrm{\times }{\mathrm{10}}^{-4}$, respectively) after introducing both Ceres and Vesta. As reported in (Laskar et al. 2011), the mutual interactions between Ceres and Vesta are so strong that the bodies cannot be considered individually when computing their orbits as well as their specific angular momentum.

Therefore, we believe this might explain the increase in the temporal variations of the longitude of the ascending node when considering only Ceres and neglecting Vesta in our computations.

After introducing the three bodies Ceres, Vesta, and Pallas together, the peak-to-peak temporal variations in Ω over the time interval [1950.0, 2050.0] continues to decrease significantly with values of $\mathrm{0}\stackrel{\mathrm{″}}{.}\mathrm{94}\mathrm{\times }{\mathrm{10}}^{-4}$ for DE405 and, $\mathrm{2}\stackrel{\mathrm{″}}{.}\mathrm{28}\mathrm{\times }{\mathrm{10}}^{-4}$ for INPOP10a.

The improvement achieved by the introduction of the three additional bodies can be easily observed in Table 10, Figs. 8(a), b for i and Figs. 9(a), b for Ω. In particular, the linear drift decreases significantly. Moreover, we note that there is a systematic oscillation with a period of 18.6 y, which should correspond to the precession period of the lunar ascending node. This residual signal is probably a consequence of the use of the EMB instead of taking into account the Earth and the Moon separately in the ephemerides.

We note that despite this oscillation, the peak-to-peak variations, with the values given above, are smaller than 1 mas for Δi and 3 mas for ΔΩ. This proves the remarkable accuracy of determination of the invariable plane.

3.2.3. Analysis of the individual contributions

We continue our analysis by evaluating the individual contributions (in terms of percentage) of each of the bodies involved in the computation of the norm of the total angular momentum. For each ephemeris and each body, we give the minimum and maximum values of this contribution, evaluated over the time interval [1950.0, 2050.0]. Results are presented in Table 8 for each of the 13 bodies of the complete system.

The contribution of Jupiter with mass M_⊙/1047.3486 and M_⊙/1047.348644, respectively, for DE405 and INPOP10a in Table 2, varies between 61.368% and 61.515% with an agreement up to the third decimal between the two ephemerides used. For Saturn of mass M_⊙/3497.888 (M_⊙/3497.9018, respectively), for the DE405 (respectively INPOP10a), the contribution varies between 24.925% and 24.957% with agreement to the fourth decimal.

For the other planets, the masses used in DE405 and INPOP10a are close to each other up to the seventh digit as can be seen in Table 2. As a consequence, we find good agreement in their respective contributions to the norm of the total angular momentum (see Table 10). As for the EMB its contribution varies between 8.5109 × 10^-2% and 8.6868 × 10^-2%. The Sun’s contribution varies between 3.4139 × 10^-4% and 1.3907 × 10^-1% with a good agreement between INPOP10a and DE405, up to the fifth decimal for the EMB, and the fourth decimal for the Sun.

We also investigated the contributions of the other bodies of the system, i.e. the two dwarf planets (134340) Pluto and (1) Ceres, as well as the two asteroids (2) Pallas and (4) Vesta. For these bodies, the masses used by the INPOP10a ephemeris are greater than those used by the DE405 ephemeris: 1% for (134430) Pluto, 1.2% for (1) Ceres, 11.3% for (2) Pallas, and 2.4% for (4) Vesta. Ceres’s contribution is about 2 × 10^-5%, and Vesta’s 5 × 10^-6% for both ephemerides, which are roughly 130 and 520 times smaller than Mercury’s contribution, respectively.

3.3. Inclinations and longitudes of the ascending nodes of the planets with respect to the invariable plane

As the fruit of a long tradition dating back to the 18th century, the ecliptic is nearly systematically chosen as the reference plane to determine the orbital elements i, Ω, $\begin{array}{}{\mathrm{˜}}_{}\\ \mathit{\omega }\end{array}$, and λ of a planet or any other moving object in the solar system (dwarf planet, asteroid, comet). Nevertheless, as we have pointed out in previous sections, the choice of adopting the ecliptic is not motivated by any physical reason, when globally studying the dynamical evolution of the solar system: it only concerns the Earth’s motion (or the EMB’s) around the Sun. We have seen in Sect. 3.2.3 that the EMB’s contribution to the total angular momentum of the solar system oscillates between 0.085% and 0.087%. Moreover, the ecliptic is a slightly moving plane and using it to define orbital parameters suggests that we choose this plane at a given epoch, i.e. J2000.0 as has it has commomly been used for approximately three decades.

In contrast, the invariable plane has two fundamental advantages inherent to its definition: first, it is fixed with respect to an inertial reference frame, as the ICRF in its more recent ICRF2 version (Ma et al. 2009). Thus, we have shown in the previous sections that its orientation with respect to the ICRF is very accurate, at least over the considered time interval [1950,2050]. Second, from its proper definition, the invariable plane represents dynamically the ideal reference plane to which one should refer to, for the orbit positioning of moving objects. In particular, it corresponds a priori to the plane for which statistically the orbital planes of moving objects are the closest, even when taking into account long-term resonant problems which should lead to high inclinations. To illustrate this property, we show in Table 8 the inclinations of the 12 bodies (planets, (134340) Pluto, (1) Ceres, (2) Pallas, and (4) Vesta) taken into account in our study when using either DE405 or INPOP10a. We note that for all the bodies (except for the Earth, of course), the inclination with respect to the invariable plane is smaller than the inclination with respect to the ecliptic. This is in particular the case for Jupiter and Saturn, for which the inclinations are $\mathrm{0}\stackrel{\mathrm{°}}{.}\mathrm{3219}$ and $\mathrm{0}\stackrel{\mathrm{°}}{.}\mathrm{9254}$ instead of $\mathrm{1}\stackrel{\mathrm{°}}{.}\mathrm{3042}$ and $\mathrm{2}\stackrel{\mathrm{°}}{.}\mathrm{4859}$, respectively.

4. Conclusions

The invariable plane of the solar system is defined as the plane perpendicular to its total angular momentum vector and passing through its barycentre. Thus, it seems to be the natural reference plane for the study of solar system bodies. It has not been the subject of any detailed studies since (Burkhardt 1982). With the adoption of a new International Celestial Reference System (ICRS) and Frame (ICRF) (see Ma et al. 1998, 2009), as well as the most recent planetary ephemerides, it seems judicious to redetermine the orientation of the invariable plane.

• 1.

  We have improved the determination of the orientation of thesolar system’s invariable plane with respect to both theorigin-equator of the ICRF and the equinox-ecliptic of the epochJ2000.0, relative to previous studies(Burkhardt 1982). To achieve this, we have usedthe two different long-term numerical ephemerides DE405(Standish 1998) and INPOP10a (Fiengaet al. 2010). Results obtained inboth cases show quite a good agreement between the twoephemerides.

• 2.

  We have studied the influence of the dwarf planets ((134340) Pluto and (1) Ceres) and the large asteroids ((2) Pallas and (4) Vesta) on the orientation and the accuracy of the determination of the invariable plane, showing a clear improvement when they are taken into account.

• 3.

  We have also computed the individual contributions of the 13 bodies considered, namely the Sun, the eight planets (the Earth being replaced by the Earth-Moon barycentre), the two dwarf planets (134340) Pluto and (1) Ceres as well the two asteroids (2) Pallas and (4) Vesta, when evaluating the total angular momentum.

• 4.

  Moreover, we have computed the inclination and the longitude of the ascending node of the planetary orbits at the epoch J2000.0, with respect to both the invariable plane and the mean ecliptic-equinox of the epoch J2000.0. In particular, we have shown that all the inclinations are smaller in the first case, except of course for the Earth.

• 5.

  Finally, our most accurate estimate of the orientation of the invariable plane with respect to the ICRF is given by an inclination of ${\mathrm{23}}^{\mathrm{◦}}{\mathrm{0}}^{\mathrm{\prime }}\mathrm{31}\stackrel{\mathrm{″}}{.}\mathrm{9}$ and a longitude of the ascending node of ${\mathrm{3}}^{\mathrm{◦}}{\mathrm{51}}^{\mathrm{\prime }}\mathrm{9}\stackrel{\mathrm{″}}{.}\mathrm{4}$. In contrast, with respect to the ecliptic-equinox of J2000.0, the orientation is given by an inclination of ${\mathrm{1}}^{\mathrm{◦}}{\mathrm{34}}^{\mathrm{\prime }}\mathrm{43}\stackrel{\mathrm{″}}{.}\mathrm{3}$ and a longitude of the ascending node of 107°34′56′′.

We assert that this determination of the invariable plane is of fundamental interest to solar system studies, in particular in the improvement of the determination of planet’s and satellites’ rotational elements as suggested by the IAU/IAG Working Group

on Cartographic Coordinates and Rotational Elements in 2006 (Seidelmann et al. 2007) and in the report of WGCCRE 2009 (Archinal et al. 2011).

---

Acknowledgments

We would like to thank Gerard Francou (SYRTE) as well as Jean Louis Simon (IMCCE) for providing us with numerical ephemeris data. We also thank Pierre Teyssandier (SYRTE) and Sebastien Bouquillon (SYRTE) for their valuable comments.

References

All Tables

All Figures