Get 20M+ Full-Text Papers For Less Than $1.50/day. Start a 14-Day Trial for You or Your Team.

Learn More →

Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework: An application to the Cassini gravitational wave experiment

Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic... Astrodynamics https://doi.org/10.1007/s42064-023-0160-x Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework: An application to the Cassini gravitational wave experiment 1,2 3 Joseph O’Leary (B), Jean-Pierre Barriot 1. School of Physics, University of Melbourne, Parkville, VIC 3010, Australia 2. Australian Research Council Centre of Excellence for Gravitational Wave Discovery (OzGrav), Parkville, VIC 3010, Australia 3. Observatoire Geo´ esique d´ de Tahiti, University of French Polynesia, BP 6570, 98702, Faa’a, Tahiti, French Polynesia ABSTRACT KEYWORDS Einstein’s theory of general relativity is playing an increasingly important role in fields general relativity such as interplanetary navigation, astrometry, and metrology. Modern spacecraft and planetary orbitography interplanetary probe prediction and estimation platforms employ a perturbed Newtonian orbit prediction and framework, supplemented with the Einstein–Infeld–Hoffmann n-body equations of motion. determination While time in Newtonian mechanics is formally universal, the accuracy of modern radiometric tracking systems necessitate linear corrections via increasingly complex and error-prone post-Newtonian techniques—to account for light deflection due to the solar system bodies. With flagship projects such as the ESA/JAXA BepiColombo mission now operating at unprecedented levels of accuracy, we believe the standard corrected Newtonian paradigm is approaching its limits in terms of complexity. In this paper, we employ a novel prototype software, General Relativistic Accelerometer-based Propagation Environment, to reconstruct the Cassini cruise-phase trajectory during its first gravitational Research Article wave experiment in a fully relativistic framework. The results presented herein agree with Received: 30 November 2022 post-processed trajectory information obtained from NASA’s SPICE kernels at the order Accepted: 27 February 2023 of centimetres. © The Author(s) 2023 1 Introduction tests of general relativity and facilitated the accurate determination of post-Newtonian parameters [3]. Of The accuracy of interplanetary navigation, planetary particular interest for weak-field tests of gravitation is the sciences, and solar system tests of general relativity has determination of the parameter γ (equal to one in general greatly benefited from modern improvements in multi-link relativity), which measures the amount of space curvature telecommunication configurations and the use of high- produced by a massive body [3]. The deflection and time frequency microwave signals between ground stations delay of photons owing to spacetime curvature can be and interplanetary probes. Precise radiometric tracking expressed in terms of γ , which is currently constrained data obtained during planetary flybys and dedicated −5 to within (2.1 ± 2.3) × 10 of the predicted unity body-centric missions have provided ample opportunity value [4]. By exploiting the multi-frequency transmission for scientific exploration of the mass, gravitational field, available at both ground station and spacecraft, the atmosphere, and orbits of celestial bodies; see Ref. [1] solar conjunction experiment performed by the Cassini and references therein for an overview of notable space probe in 2002 [5] provided the first complete solar plasma missions dedicated to planetary sciences. calibration for both uplink and donwlink signals [6]. The The use of joint Ka/X band microwave links coupled with next-generation multi-frequency technologies [2] Cassini radioscience experiment improved the estimation has yielded unprecedented accuracy for solar system of γ by approximately two orders of magnitude, which B joe.oleary@unimelb.edu.au 2 J. O’Leary, J.P. Barriot −3 was previously constrained to within ∼ 10 of the general [23] correctly account for the relative motion of the relativistic prediction using data from the Viking Martian Sun with respect to the solar system’s barycenter. The lander relativity experiment [7]. expected theoretical difference affects the obtained value −5 Through a series of resolutions, the International of the post-Newtonian parameter γ = (2.1 ± 2.3) × 10 −4 Astronomical Union (IAU) [8] suggests that all problems by the amount ± 1.2 × 10 owing to the barycentric in the field of astronomy or astrodynamics should be velocity of the Sun [2, 21]. Nevertheless, the legacy formulated within the framework of Einstein’s theory software used in the Cassini radioscience experiment, of general relativity. These include: (i) the derivation i.e., Orbit Determination Program [13, 14], employed an of dynamical equations of motion; (ii) the definition approximate model for the motion of the Sun, with a of observables; and (iii) the transformations between negligible effect on the final result [20]. reference frames and time scales. The high precision It is expected that the current best estimate range and Doppler observables combined with data from of γ will be improved in the near future by the onboard accelerometers [9–12] necessitate that modern Mercury Orbiter Radioscience Experiment (MORE) orbitography software employs high-fidelity dynamical onboard the ESA/JAXA BepiColombo mission [24]. and measurement models for spacecraft propagation and The radioscience instrumentation onboard BepiColombo nonlinear state and parameter estimation, respectively will provide plasma-free range information in addition [13, 14]. to stable Ka/X band microwave links [25], while Popular orbitography software such as the French simultaneously benefiting from the direct measurement Space Agency Geo´ desie ´ par Inegrations t´ Numeriques ´ of non-gravitational forces in the local frame of the Simultanees ´ (GINS) [ 15] or NASA’s Mission Analysis, spacecraft by an onboard accelerometer [26]. The Operations, and Navigation Toolkit Environment novel 24 Mcps pseudo-noise ranging instrumentation (MONTE) [16] currently describe the motion of onboard BepiColombo [27] will provide two-way range interplanetary spacecraft via a classical Newtonian information at the centimeter level, so that second- framework linearly corrected with effective forces, order post-Newtonian contributions in the light-time accounting for the effects of general relativity via the so- solution need to be considered in the data reduction called n-body Einstein–Infeld–Hoffmann (EIH) equations process [28]. In addition, future space-based gravitational of motion [17]. Given the stringent accuracy requirements wave detectors, e.g., the Laser Interferometer Space associated with fields such as astrometry, geodesy, Antenna (LISA), are expected to detect the propagation −5 celestial mechanics, and metrological sciences [18], of low-frequency gravitational waves (10 –1 Hz), spacecraft and interplanetary probe orbit propagation complementing the scientific operations of ground- tools based on the so-called “Newton + correction” based detectors such as the Laser Interferometer framework need to include subtle relativistic effects Gravitational wave Observatory (LIGO), VIRGO, and in both state and measurement models to reflect KAGRA collaborations, which focus on higher frequencies the rapid improvements in modern measurement in the gravitational wave frequency spectrum. LISA technologies. The aforementioned ad-hoc process of employs time-delay interferometry (TDI) to reduce reducing relativistic phenomena to simple perturbations the effects of laser frequency and optical bench noise or relativistic corrections introduces a neo-Newtonian [29]. TDI analysis employs solar system barycentric [19] interpretation of the theory of general relativity, coordinate time while each LISA spacecraft records data which can lead to spurious results when equivalence according to the individual spacecraft proper time [30]. between Newtonian and relativistic concepts is assumed. Hence, sophisticated four-dimensional general relativistic Moreover, the “Newton + correction” approach creates spacetime transformations are required for TDI data the likelihood of neglecting potentially important reduction to detect gravitational waves using space-based relativistic contributions if care is not taken. Indeed, it interferometry. Bearing this in mind, we believe that has been discussed at length in relevant literature [20–22] the classical “Newton + correction” paradigm is now whether the dynamical equations and the corresponding reaching its limits in terms of complexity. measurement models used in the data reduction process In a recent article [31], we presented the first results of of the celebrated Cassini solar-conjunction experiments a prototype software, namely General Relativistic Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 3 Accelerometer-based Propagation Environment 2 GRAPE: A brief overview (GRAPE), which describes the perturbed motion In this section, we discuss the mathematical framework of Parker Solar Probe [32], Mercury Planetary Orbiter implemented in GRAPE. First, we provide a high-level [33], and Molniya-like [34] interplanetary probes and overview of the platform and discuss some aspects of spacecraft subject to a radiation-like non-gravitational the software that are unique to general relativity and four-force using the complete framework of general not common with classical propagation environments. relativity. GRAPE employs extended relativistic In Section 2.1, we present the equations of motion and equations of motion that account for non-gravitational discuss the treatment of non-gravitational perturbations forces using either end user supplied accelerometer data in a fully relativistic framework. or approximate dynamical models in the local frame In the framework of general relativity, the unperturbed of the spacecraft. We exploit the tetrad formalism [35] motion of test particles is described using a geodesic to compare dynamical quantities and non-gravitational equation parameterized by an arbitrary affine parameter forces in both co-moving and global frames. The s. For massive test particles such as interplanetary probes novel approach adopted by GRAPE assumes that and near-Earth spacecraft, it is common to parameterize the local structure of spacetime is unperturbed by equations of motion using their associated proper time τ external non-gravitational forces, and that all associated such that s ≡ τ . While proper time parameterization is gravitational and relativistic contributions are contained convenient when working in a theoretical framework, it is and accounted for at the order of the spacetime metric inconvenient for practical purposes. Spacecraft operators, tensor. In this paper, we discuss the preliminary steps tracking stations, and orbit determination as well as data taken to extend GRAPE to the operational planning reduction software reference the state, error covariance, of future missions within the solar system. Further, and observables with respect to a derived external we describe the initial approach adopted to interface coordinate time scale t such as terrestrial, barycentric, it with NASA’s SPICE kernels [36], and verify our or ephemeris times. Hence, GRAPE internally adopts results by comparing the pre-fit Doppler residuals 0 a global coordinate time t = x /c during numerical using GRAPE and SPICE with high-precision Ka/Ka propagation, and the corresponding equations of motion band Doppler data obtained from the first Cassini are transformed from an affine to a coordinate time probe gravitational wave experiment (GWE1) [5]. parameterization accordingly. The remainder of this paper is organized as follows. The incorporation of non-gravitational forces in In Sections 2 and 3, we summarize the equations of classical Newtonian propagation environments is a motion and associated numerical integration procedures well-recognized problem and remains an active area implemented in GRAPE, respectively. In Section 4, we of research in operational astrodynamics [37]. As present new results that extend GRAPE for practical discussed at length in Ref. [31], the treatment of mission planning using radiometric observables and non-gravitational perturbations in general relativity determine that GRAPE recovers the Cassini trajectory has not received the same attention in the literature at the order of centimetres. We conclude the paper in and requires further elucidation. The nonlinear field Section 5, providing a brief discussion on future plans equations of general relativity [38, 39] can be decomposed for the platform. according to contributions from the Einstein tensor Unless otherwise stated, we assume the Einstein G and the energy momentum tensor T . The µν µν summation convention for repeated indices; the use Einstein tensor is composed of the metric tensor of Latin indices (i, j = 1, 2, 3) is reserved for spatial and quantities derived from it, e.g., the Ricci and coordinates, and Greek indices (α, β = 0, 1, 2, 3) denote Riemann tensors, which encode the geometry of spacetime spacetime coordinates, with the 0th component reserved so that all gravitational phenomena are captured α i for time with x = (ct, x ). Finally, we assume a entirely at the order of the spacetime metric tensor. spacetime metric signature according to (+,−, −, − ), Similarly, quantities derived from the energy-momentum a˙ denotes the derivative of a with respect to time, and c tensor encode information related to non-gravitational denotes the speed of light in vacuum. phenomena. In the parlance of operational astrodynamics 4 J. O’Leary, J.P. Barriot i i and interplanetary navigation, gravitational and non- where we introduced the dimensionless quantity β = v /c i i gravitational perturbations are obtained from the metric in terms of the test particle four-velocity v = dx /dt, and and energy-momentum tensors, respectively. Hence, the the Christoffel symbols of the second kind are denoted perturbed equations of motion employed by GRAPE are by Γ . It is worth noting that the components of the βγ derived from the divergence of the energy-momentum four-force in Eq. (5) are given by Eq. (3), and in the tensor, which is linearly decomposed according to appropriate Newtonian limit, Eq. (5) reduces to Newton’s contributions from the energy-momentum of an isolated second law perturbed by an external force f , which can test particle and the corresponding stresses acting on be most easily observed via the explicit calculation of the the system that produce geodesic and perturbed motion, Christoffel symbols given by Eq. (2) in Ref. [42]. On close inspection, Eq. (3) raises two important respectively [40, 41]. practical questions: (i) How are the components of K 2.1 GRAPE: Equations of motion obtained? (ii) Given that non-gravitational forces are To account for non-gravitational forces in Einstein’s measured in the local frame of the spacecraft, how does theory of general relativity, we decompose the energy- GRAPE relate the quantities between local and global momentum conservation equation according to frames, such as the IAU-recommended [8] geocentric and µν µν µν barycentric celestial reference frames? We address both ∇ T ≡ ∇ (Θ + S ) = 0 (1) µ µ questions below. µν where we introduced the stress tensor S , and we The inner product of two four-vectors is a scalar µν interpret Θ as the energy-momentum of an isolated invariant I. In particular, it can be readily shown that µν test particle. In its simplest form, Θ is given in terms in all reference frames, µ µ of mass density ρ and four velocity u ≡ dx /dτ , with µ ν dx dx µν µ ν I ≡ g = c (6) Θ ≡ ρu u . Equation (1) gives rise to two distinct µν dτ dτ µν µν cases: (i) S = 0 and (ii) S ̸= 0. In the absence Differentiating Eq. (6) with respect to proper time and of stresses, it is well established that Eq. (1) yields a substituting the expression given by Eq. (3), we find that µν geodesic equation of motion. In the case where S = ̸ 0, at every point in spacetime, the four-force is orthogonal we find to the four-velocity with µν ν ∇ Θ = ρK (2) µ σ 2 µ f u ≡ K − (K u ) u /c u = 0 (7) µ µ σ µ µν ν where ∇ S = −ρK is the force density associated For Eq. (7) to be satisfied in the local (L) co-moving with an arbitrary non-gravitational four-force f . µ frame of the spacecraft where u = (c, 0, 0, 0), we find Evaluating the term on the left-hand side of Eq. (2) α i that f = 0, K . Hence, in the local co-moving frame, yields [31]: the non-gravitational perturbations are purely spatial and µ ν µ µν µ ν 2 f ≡ u ∇ u = K g − u u /c (3) ν ν are provided to GRAPE by end users using either known dynamical models or data obtained from an onboard where it is clear that geodesic motion is recovered accelerometer. The three-dimensional (Newtonian) non- when K ≡ 0. As previously discussed, Eq. (3) can gravitational forces acting on interplanetary probes be equivalently expressed with respect to a non-affine used in classical orbitography are well known and coordinate time t following several applications of the documented; see Ref. [43] for explicit examples. However, chain rule: −1 α α in the framework of general relativity, the corresponding dx dx dt models for the impulsive and continuous four-forces dt dτ dτ (4) −2 2 α 2 α 2 α are not readily available, posing a theoretical obstacle d x dt d x d t dx = − 2 2 2 for general relativistic orbitography. By modeling the dt dτ dτ dτ dt non-gravitational four-forces in the local frame of the so that the final equations of motion employed by GRAPE spacecraft first, we circumvent this issue, whereupon the are given by 2 i global four-forces can be obtained using the so-called 1 d x i i k i j k 0 i 0 i k = −Γ − 2Γ β − Γ β β + Γ β + 2Γ β β 00 0k jk 00 0k 2 2 tetrad formalism [35]. c dt To compare dynamical quantities between local and 1 dτ 0 i j k 0 i i + Γ β β β − f β − f (5) jk 2 global frames, GRAPE employs the tetrad formalism c dt Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 5 as follows. At every point in spacetime, we construct a admits first integrals or quadratic invariant quantities, pseudo-Cartesian reference system by introducing e symplectic integration procedures ensure that these (ν ) which satisfy remain constant (or bounded in some cases) throughout ν (α) (β ) the numerical integration. δ = g e e , g = e e η (8) αβ µν µν αβ µ ν (α) (β ) The dynamical systems associated with classical We note that e consists of four orthonormal vectors as (ν ) orbitography software are both time- and velocity- opposed to a single tensor quantity known collectively as dependent (e.g., atmospheric drag perturbations), so a tetrad or vierbein [44]. The vector of interest is indicated that the resulting system of equations are not readily by the subscript in parentheses and the component of integrable. Hence, in a classical framework, structure- that vector is denoted by the superscript. preserving integration schemes do not provide an An important feature of the tetrad formalism is advantage over non-symplectic algorithms. This is not that quantities projected onto the tetrad have physical the case for general relativistic orbitography. Given that meaning. Therefore, Eq. (6) is valid in every reference frame, regardless of (0) (0) µ (i) (i) ν dx = e dx , dx = e dx (9) µ ν time- or velocity-dependent external non-gravitational denote locally measured spatial and temporal intervals four-forces, a unique opportunity to exploit structure- [35]. From Eq. (9), it is clear that the force components preserving integration schemes is available with general are obtained according to relativistic mission planning platforms. GRAPE employs (σ ) (σ ) µ ν ν (σ ) f = e f , f = e f (10) the entire suite of numerical integration procedures µ (σ ) presented in Ref. [46], where the reader is referred Equation (10) provides a direct transformation between to Refs. [31, 46, 47] for the specific Butcher tableau local and global force components when e denotes the (ν ) coefficients required for practical implementation. In components of the co-moving tetrad. If, for example, this section, we review implicit Runge–Kutta algorithms e denoted the components of the so-called natural (ν ) and present a derivation of their structure-preserving tetrad, then a further Lorentz transformation is required properties for quadratic invariants. Finally, we note that to account for the relative motion (see Ref. [31, Section this section is based on Ref. [31, Section 4] and remains 2.4] for further details). largely unchanged from its original presentation. 3 GRAPE: Symplectic integration 3.1 GRAPE: Implicit Runge–Kutta proce- scheme dures In general, classical numerical integration algorithms do not preserve the first integrals (constants of Implicit Runge–Kutta algorithms proceed as follows. motion) associated with dynamical systems; thus, over Given an arbitrary system of N differential equations long integration timescales, energy dissipation occurs, with associated initial conditions: leading to incorrect qualitative dynamical behavior. dx i j i i = f (t, x ), x (t ) = x (i, j = 1, 2, . . . , N) Symplectic integration schemes [45] are advanced dt (11) numerical procedures that preserve the symplectic the solution is propagated x → x with an integration n n+1 structure of Hamiltonian systems such that the so-called step-size h according to flow of the system defines a canonical transformation and belongs to a class of numerical methods known i i i x = x + h b k (12) n+1 n j as structure-preserving geometric integrators. In the j=1 celebrated work of Butcher [46], the author introduces i i i i and discusses the properties of implicit Runge–Kutta k = f t + c h, x + h a k (13) n j jl j n l l=1 schemes, presenting the Butcher tableau coefficients for a where Eq. (13) is defined implicitly. Presently, GRAPE series of s-stage, nth-order algorithms, known collectively determines the individual k using fixed-point iteration, as Gauss collocation methods [45]. Implicit Runge–Kutta schemes are symplectic by definition (see Section 3.2 for where the initial value is assumed to be zero everywhere further details). Hence, given a dynamical system that so that (k ) = 0. Improved strategies will be employed j 6 J. O’Leary, J.P. Barriot in future work to initialize k , which requires further where Eq. (20) is obtained by squaring Eq. (12). investigation [45]. Furthermore, x denotes the ith component of an arbitrary vector at time t = t , and k denotes the mth 3.2 Quadratic invariants approximate slope associated with the ith differential equation. For Eq. (16) to be classified as a first integral, General relativity gives rise to a unique dynamical we must have Q = Q = . . . = Q = C,∀n, where quantity (see Eq. (6)) formally known as a quadratic n+1 n 0 invariant [45]. Thus, the symplectic integration Q ≡ Q(x ) denotes the value of the quadratic function procedures implemented in GRAPE provide an additional (16) at time t and C is an arbitrary constant. Hence, layer of verification during the propagation interval. The using Eq. (13), we can equivalently express Eq. (20) coefficients b and a identified in Ref. [ 46] satisfy the according to i ij constraint: i i j i j M x x = M x x + 2h b M X f ij ij ℓ ij n+1 n+1 n n ℓ b a + b a − b b = 0 (i, j = 1, 2, . . . , s) (14) i ij j ji i j ℓ=1 s s X X which consequently indicates that Gauss collocation 2 i j + h (b b − b a − b a ) M k k (21) ℓ m ℓ ℓm m mℓ ij ℓ m methods [46] are symplectic and completely preserve ℓ=1 m=1 quadratic invariants. The proof is outlined below, and i i i where, for brevity, we have introduced k = f (T, X ), the reader is referred to Ref. [45, Chapter 4] for the s i i i with T = t + c h and X = x + h a k . n j jl n n l=1 l complete details. According to Eqs. (17) and (14), the second and final A first integral I(x ) associated with an arbitrary terms on the right-hand side of Eq. (21) are zero; system of differential equations (11) satisfies the condition therefore, for Runge–Kutta methods satisfying Eq. (14), in Eq. (15): j i i j M x x = M x x . Thus, with respect to the ij ij n+1 n+1 n n dI ∂I α α α i j current problem, we find Q(x ) = I(η (x ), u ), where µν = f (t, x ) = 0 (15) dt ∂x I is defined according to Eq. (6), and the condition given We define a quadratic function of Eq. (11) according to by Eq. (19) corresponds to the orthogonality condition k i j Q(x ) = M x x (16) ij associated with the four-velocity and four acceleration; see Eq. (7). It is important to note that in order to enforce where M is an arbitrary symmetric matrix. Hence, ij strict symplecticity, we must express I in the local co- Q(x ) is the first integral of Eq. (11) if Eq. (17) holds moving tetrad so that at each numerical integration step, [45]: the spacetime is assumed to be flat, which is a formal i j k M x f (t, x ) = 0 (17) ij requirement according to the condition dM /dt = 0. ij where it is clear that we must have dM /dt = 0 which ij To conclude this section, we reiterate the important is demonstrated as Eq. (18): point: In classical orbit problems, symplectic integration i j d dM dx dx ij procedures have been shown to have attractive numerical i j i j j i M x x = x x + M x + M x ij ij ij dt dt dt dt properties for the long-term numerical integration of = 0 (18) Hamiltonian systems [45]. Implicit Runge–Kutta schemes and using the symmetry M = M we find [46] are not restricted to Hamiltonian systems. They ij ji are a family of implicitly defined numerical procedures dx i i j M x = M x f = 0 (19) ij ij for the numerical solution of any system of differential dt if and only if M is a constant. equations, e.g., stiff and non-stiff problems [ 45]. The ij According to Eq. (12), we obtain Eq. (20): Gauss collocation methods introduced above are a particular example of an implicit Runge–Kutta scheme j j i i j i M x x = M x x + h b M x k ij ij ℓ ij n+1 n+1 n n n with the added property that their Butcher tableau ℓ=1 coefficients satisfy the constraint given by Eq. (14). i j + h b M k x Hence, they naturally preserve quadratic invariants, as m ij m n m=1 discussed above. s s X X In many classical problems, the introduction of 2 i j + h b b M k k (20) ℓ m ij ℓ m velocity- and time-dependent perturbations results in ℓ=1 m=1 Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 7 ④ ⑤ a system of differential equations that may not readily Doppler ; and (iii) using the SPICE spkezr program, admit first integrals, such that the structure-preserving we perform a direct comparison of the GRAPE ephemeris properties of symplectic integration schemes can no longer data with the Cassini SPICE kernel and project the be exploited. The equivalent problem in the framework residuals onto the radial, along-, and cross-track axes in of general relativity is different; regardless of the specific the local time-dependent spacecraft coordinate system; non-gravitational four-forces entering Eq. (5), the inner see Section 3.3.3 in Ref. [37] for further details. In the product of two four vectors will always be an invariant. In present analysis, the non-gravitational perturbations particular, Eq. (6) is always satisfied, despite the motion entering Eq. (23) are inferred directly from the Cassini of the spacecraft or probe being subject to external forces. SPICE kernel, and we assume that they are held constant throughout the experiment. Given that this is the first Accordingly, we choose to employ an implicit Runge– application of GRAPE to real radiometric tracking data, Kutta scheme for the present analysis, which ensures we persevere with this assumption and note that it is that Eq. (6) does not diverge from c beyond machine an approximation. We further note that the GRAPE precision, providing an additional accuracy check on estimation module is currently being developed, where modifications to the platform. While the results presented specific time-dependent non-gravitational force models in Section 4 are restricted by the duration of GWE1, i.e., will be applied to parameter estimation problems for O (days), we note that GRAPE can be used for long-term probes during interplanetary cruise. The explicit details integration, e.g., generating solar system ephemerides and steps required to reproduce the present analysis are [O (centuries)] [48, 49], which is of interest to the pulsar discussed below. Finally, we note that by construction, timing community for the detection of gravitational the routines utilized internally by the SPICE toolkit waves using time of arrival residuals [50]. The necessary automatically take into account the Earth ephemeris, modifications to GRAPE are minor and will be explored the location of Earth-based tracking stations, and the in a future study. Earth’s rotation. For approximately 40 days starting November 26, 4 Results 2001, the Cassini probe was tracked during solar In this section, we demonstrate the modifications made opposition in search of gravitational waves in the to GRAPE since Ref. [31] in preparation for future millihertz frequency band [52]. Given that gravitational practical mission planning, which we verify through three wave detection methods depend on the timescale of independent strategies: (i) using GRAPE, we generate the radiation, ground-based detectors based on laser an ephemeris for the Cassini probe for the duration of interferometry cannot separate signals from seismic −5 GWE1 and ensure that Eq. (6) is preserved throughout and other noise sources in the frequency range 10 – −2 the integration period; (ii) using the SPICE mkspk 10 Hz. Consequently, searches for low-frequency program, we construct a SPICE kernel containing the gravitational waves must employ space-based detectors GRAPE-generated Cassini trajectory. We calculate the [5]. The gravitational wave experiments performed light-time solution using the GRAPE and Cassini SPICE by Cassini were, at the time, the most sensitive ②③ kernels at each observation epoch and compare the Doppler tracking campaigns achieved using interplanetary corresponding values with the observed GWE1 two-way probes. Moreover, at a distance of approximately ∼ 8 AU, Cassini became the largest gravitational wave ① See https://naif.jpl.nasa.gov/pub/naif/toolkit docs/FORTRA N/ug/mkspk.html for further details. antenna ever used [5]. The success of the Cassini ② See https://naif.jpl.nasa.gov/pub/naif/pds/data/co-s j e v-spi gravitational wave experiments and solar system tests of ce-6-v1.0/cosp 1000/data/spk/ to download the relevant Cassini general relativity was mediated by two major technical SPICE kernel and https://naif.jpl.nasa.gov/pub/naif/pds/data/ co-s j e v-spice-6-v1.0/cosp 1000/data/spk/041014R SCPSE 0 upgrades. First, Cassini was equipped with a novel 1066 04199.lbl for specific kernel details. multi-link telecommunication configuration with two ③ It is important to note that the trajectory information contained in individual Cassini SPICE kernels are not observables. Rather, they are the result of an orbit determination fit using two-way ④ The dataset ID used here is CO-X-RSS-1-GWE1-V1.0. It can be Doppler and range observables obtained from Cassini’s two low downloaded at NASA’s Planetary Data System. gain antennae and high gain antenna; see Ref. [51] for explicit ⑤ See https://naif.jpl.nasa.gov/pub/naif/toolkit docs/FORTRA details on reconstructing the Cassini trajectory. N/spicelib/spkezr.html for further details. 8 J. O’Leary, J.P. Barriot 2 2 2β − 1 µ s˙ s˙ uplink carriers in the X and Ka bands as well as three k i j − + γ + (1 + γ ) c r c c downlink carriers in the X, Ka1 (coherent with the X- jk k̸=j uplink), and Ka2 (coherent with the Ka-uplink) bands, i j j j 2 (1 + γ ) 3 x x˙ − x x˙ i j − x˙ x˙ − considerably reducing the electromagnetic wave phase 2 2 c 2c r ij scintillation owing to solar plasma [6, 52]. Second, the j j i j radioscience experiments performed by Cassini benefited + x x¨ − x x¨ 2c from a dedicated advanced media calibration system 1 µ i j i + x − x (2 + 2γ ) x˙ located at the Goldstone Deep Space Communications 2 3 c r ij j̸=i Complex (DSS-25). The novel system consisted of water 3 + 4γ µ x¨ vapor radiometers, pressure sensors, and microwave j i j − (1 + 2γ ) x˙ x˙ − x˙ + 2c r ij temperature profilers that calibrated the frequency j̸=i modulation caused by the wet and dry components of 1 dτ 0 i i − f β − f (23) the troposphere [5, 6, 52]. c dt The curved spacetime associated with the celestial i i where we introduced the notation x x = x · x to −5 bodies is described (up to order c in the time denote the inner product between two vector quantities, −3 component, up to order c in the spatial components, and we assume that the components of x are given −4 and up to order c in the off-diagonal components) by in terms of Cartesian coordinates. Equation (23) is the components of the n-body metric tensor (see Ref. [14, defined implicitly in terms of x¨ , such that the formal Section 2]): numerical solution is obtained by iteration. However,   X X µ x˙ for the present demonstration, we assume that terms 2 + 2γ j 2γ µ j j   g = , g = −δ 1 + , 0k kl kl i 2 3 2 c r c r of the form x¨ /c can be replaced by their Newtonian ij ij j̸=i j̸=i counterparts. It is important to note that the final term   X X X µ s˙ 2 µ 2β µ 1 + 2γ appearing in Eq. (23) completely captures the general j j j   g = 1− + − 2 4 4 c r c r c r relativistic contributions of an arbitrary non-gravitational ij ij ij j̸=i j̸=i j̸=i four-force f . X X X 2 (2β − 1) µ µ 1 ∂ r j k ij + − µ The calculation of observable quantities in a complete 4 4 2 c r r c ∂t ij jk j̸=i k̸=j j̸=i general relativistic framework [53] relies on highly (22) accurate light-time solutions and is profoundly different where the gravitational parameter of the jth celestial from classical signal-processing paradigms [54, Chapter body is given by µ ; the kth component of its coordinate 2]. In the present study, light-time is calculated using the k 2 k velocity is given by x˙ ; s˙ = x˙ ; β and γ denote j j j approach presented in Ref. [14, Chapter 8]. We note that post-Newtonian parameters; the relative barycentric light-time solutions obtained by numerically solving the distance between bodies i and j is given by r ; and we ij null geodesic equations where the effects of solar plasma have introduced the Kronecker delta function according are contained within a modified metric tensor [ 55] are to δ . In the present demonstration, we assume that kl currently being explored. body i refers to the Cassini probe and bodies j = Figure 1 displays the temporal evolution of the 1, . . . , n − 1 (j ̸= i) refer to the Sun, the nine planets, quadratic invariant I given by Eq. (16). The quadratic and the Moon, where the initial conditions of all n bodies invariant I depends explicitly on the components of the are obtained from the appropriate SPICE kernels on Nov α numerically integrated four-velocity u and implicitly 26, 2001 at 05:00:00 Barycentric Dynamical Time using α on the numerically integrated spacetime coordinates x the SPICE spkezr program. The extended equations of through the metric tensor components g given by µν motion implemented in GRAPE follow directly from the Eq. (22). At each integration step, GRAPE internally definition of the metric tensor components and Eq. (5), checks whether I deviates from machine precision. and are given by Eq. (4-26) in Ref. [14, Section 4]: This is employed as an additional accuracy check j i 2 i X X d x µ x − x 2 (β + γ ) µ l on modifications made to GRAPE. It is clear from = 1 − 2 2 dt r c r il Fig. 1 that the numerically determined four-velocity ij j̸=i l̸=i Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 9 −16 1.0 ×10 SPICE 2.0 GRAPE 0.5 1.5 1.0 0.0 0.5 −0.5 0.0 −1.0 −0.5 −1.0 −1.5 −1.5 −2.0 −2.0 0 5 10 15 20 25 30 35 40 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 Time (day) Time (day) Fig. 1 Time evolution of the quadratic invariant for the Fig. 2 Doppler residuals using SPICE (black) and GRAPE Cassini probe during GWE1. The orbit is propagated using (cyan). The GRAPE ephemeris is propagated using an an integration time-step h = 300 s with an initial condition integration time-step h = 300 s with an initial condition obtained from NASA’s SPICE kernels on Nov 26, 2001 at obtained from NASA’s SPICE kernels on Nov 26, 2001 05:00:00 Barycentric Dynamical Time. The time axis is given at 05:00:00 Barycentric Dynamical Time. The SPICE in days since the initial epoch. We define δI according to ephemerides are pseudo-observations deduced using a least- 2 2 δI ≡ (I − c )/c , where I is determined at each numerical squares fit from radiometric observables. The time axis is integration step and is defined in Eq. (16). GRAPE preserves given in days since the initial epoch. −16 I in the order of 10 . in orbit determination problems, the goal is to obtain components satisfy Eq. (6) at machine precision, which Gaussian residuals with zero mean. In this regard, the −32 can be further reduced to 10 by selecting 128-bit results presented herein are slightly different. In its (quadruple) precision internally at a computational simplest form, orbit determination software is typically cost. The computational demand for double (64-bit decomposed into two modules: prediction and estimation. precision) and quadruple precision is O (seconds) and The goal of the present study is not estimation via least- O (minutes), respectively. We note that a thorough squares or Kalman filter analysis. Instead, this study aims analysis benchmarking the GRAPE platform is deferred to present the first results of the GRAPE platform during for future study and is outside the scope of the present the prediction phase. As noted above, the estimation paper. module of GRAPE is currently under development and Figure 2 compares the pre-fit Doppler residual using will be the topic of a future paper on general relativistic GRAPE and Cassini SPICE trajectory kernels, as state and parameter estimation. discussed in Steps (4)–(6). The details of the Cassini Figure 3 compares the GRAPE and Cassini ˆ ˆ ˆ kernel and observational data used in the present analysis ephemerides in the time-dependent (R, S, W ) spacecraft are presented in Footnotes ② and ④. For clarity, we frame, as discussed in Step (7). Using GRAPE, the plotted only the residuals for the first pass of Cassini reconstructed Cassini trajectory during cruise agrees −4 during GWE1. However, we note that the difference with the SPICE ephemeris at the level of 10 km. in the residuals remains at the ∼ 0.2 Hz level, with The projected residuals in the cross- and along-track no noticeable secular variation during the first 1 × 10 directions, respectively denoted as ∆ W and ∆ S, are observations. The pre-fit Doppler residuals are slightly negligible. The projected residual along the radial biased owing to uncertainties in the initial state of the direction, denoted as ∆ R, grows secularly over the spacecraft obtained from SPICE. The bias is corrected in duration of the Cassini GWE1. During cruise, the the estimation phase, e.g., by employing a (sequential) dominant non-gravitational force is in the radial direction. batch least-squares or Kalman filter algorithm [ 37]. The Given that we assume the non-gravitational forces remain initial-state bias in Fig. 2 is standard and is comparable to constant during GWE1, the projected residual observed that of other modern orbit prediction and determination along the ∆ R axis is expected. We conjecture that software; see Fig. 1 in Ref. [1] for example. Typically, once specific time-dependent non-gravitational four-force δI Residual (Hz) 10 J. O’Leary, J.P. Barriot at each observation time using finite differences. ∆R 0.0004 (5) Define the non-gravitational (NG) perturbations ∆W acting on the probe according to ∆S 0.0003 a ≃ a − a (24) NG S G The GRAPE (G) and SPICE (S) acceleration vectors 0.0002 are respectively defined by 0.0001 a ≃ a , a ≃ a + a (25) G M S M NG where we have introduced the acceleration owing 0.0000 to the n-body metric tensor components only as a . Under the assumption that a ≃ a during M NG 0 0 5 10 15 20 25 30 35 40 Time (day) the short duration cruise phase (< 40 days), where a denotes a constant vector, we employ Eq. (10) Fig. 3 Comparison of SPICE and GRAPE ephemerides’ residuals. The GRAPE ephemeris is propagated using an at each integration time-step and produce a new integration time-step h = 300 s with an initial condition GRAPE-generated kernel to be used with the SPICE obtained from NASA’s SPICE kernels on Nov 26, 2001 software. The ephemerides associated with the new at 05:00:00 Barycentric Dynamical Time. The SPICE kernel now take into account the non-gravitational ephemerides are pseudo-observations deduced using a least- squares fit from radiometric observables. The time axis is forces acting on the probe under the assumption that given in days since the initial epoch. they remain constant throughout the integration procedure. Again, we reiterate that the adopted models are implemented in GRAPE, the ∆ R residual ⑥ approach is an approximation. The implementation will be significantly reduced . of specific time-dependent models in the local tetrad It is important that the results presented in this is postponed to future work. manuscript are easily reproducible. As a result, and in (6) For Fig. 2, we calculate at each observation time the anticipation of a future public release of the GRAPE two-way light-time solution between the spacecraft software, we outline the specific steps required to and Earth antennae (see Ref. [14, Chapter 8] for reproduce Figs. 1–3: explicit details) and compare the calculated Doppler (1) Generate an ephemeris using GRAPE by numerically frequency with the observed value. solving Eq. (23) subject to the condition that α (7) For Fig. 3, we introduce a time-dependent local f ≡ 0. This step produces an ephemeris for the ˆ ˆ ˆ spacecraft coordinate system (R, S, W ) given in Cassini trajectory during the first gravitational wave terms of the probe position r and velocity vector v experiment under the assumption that the motion as Eq. (26) [37]: of the probe is governed entirely by the n-body solar r r × v system metric tensor components given by Eq. (22). ˆ ˆ ˆ ˆ ˆ R = , W = , S = W × R (26) r |r × v| (2) Using the SPICE program mkspk, create a SPICE We project the residual position vector between kernel from the GRAPE-generated ephemeris GRAPE and SPICE, denoted by δ r ≡ r − r , along G S obtained in Step (1). each of the axes in the local spacecraft frame so that (3) Retrieve the approximate probe coordinates at the residuals plotted in Fig. 3 are given explicitly by each observation time relative to the solar system ˆ ˆ ˆ barycenter using both GRAPE-generated and ∆R ≡ δ r · R, ∆W ≡ δ r · W , ∆S ≡ δ r · S (27) Cassini SPICE kernels, using the SPICE program spkezr. 5 Discussions and conclusions (4) Given the SPICE and GRAPE trajectory We presented the first practical results of a novel general information, numerically calculate the acceleration relativistic orbitography platform, GRAPE. While the ⑥ We performed several experiments with h varying between 15 initial results are in good agreement with NASA’s and 300 s. The results of each experiment were comparable, SPICE software, further work is required. It is of with no noticeable improvement in the solution as h decreases. Accordingly, we adopted h = 300 s. utmost importance that GRAPE can accurately provide Projected residual (km) Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 11 light-time solutions for the calculation of radiometric Research Council Centre of Excellence for Gravitational observables in a fully relativistic framework. Rather Wave Discovery (OzGrav) (Grant No. CE170100004). than adopting the standard approach in the literature, JPB was funded by a DAR grant in planetology from the where light-time is approximated by a series of linearly French Space Agency (CNES). superimposed post-Newtonian corrections that gradually Funding note increase in complexity and are error-prone for higher- order corrections, we are currently implementing light- Open Access funding enabled and organized by CAUL time solutions in GRAPE by numerically solving the null and its Member Institutions. geodesic equation given in Ref. [55]. The latter approach guarantees that all relevant relativistic contributions to Declaration of competing interest the dynamics are appropriately accounted for at the order The authors have no competing interests to declare that of the metric tensor and voids the possibility of neglecting are relevant to the content of this article. potentially important terms [20–22]. The motivation for designing a general relativistic References orbit prediction and estimation tool is not to replace [1] Lari, G., Schettino, G., Serra, D., Tommei, G. Orbit the current popular orbitography platforms based on the determination methods for interplanetary missions: corrected Newtonian approach. Rather, it is to provide Development and use of the Orbit14 software. a fully consistent and highly accurate benchmark for Experimental Astronomy, 2022, 53(1): 159–208. current and future platforms. With radioscience data [2] Bertotti, B., Comoretto, G., Iess, L. Doppler tracking preprocessing software now publicly available (see for of spacecraft with multi-frequency links. Astronomy and example Ref. [56]), archival tracking data from a myriad Astrophysics, 1993, 269(1–2): 608–616. of space missions can be easily reduced, thereby providing [3] Will, C. M. Theory and Experiment in Gravitational Physics. Cambridge University Press, 2018. ample opportunity to revisit historical missions and [4] Bertotti, B., Iess, L., Tortora, P. A test of general verify future modifications to the tool. We expect the relativity using radio links with the Cassini spacecraft. current manuscript to be the first in a series of papers Nature, 2003, 425(6956): 374–376. describing new approaches to general relativistic orbit [5] Kliore, A. J., Anderson, J. D., Armstrong, J. W., Asmar, determination (see for example Ref. [47, Chapter 7]) S. W., Hamilton, C. L., Rappaport, N. J., Wahlquist, H. applicable to missions requiring high precision, e.g., D., Ambrosini, R., Flasar, F. M., French, R. G., et al. BepiColombo [12], JUICE [57], and VERITAS [58, 59] Cassini radio science. In: The Cassini-Huygens Mission. missions, as well as the space-based gravitational wave Dordrecht: Kluwer Academic Publishers, 2005: 1–70. [6] Tortora, P., Iess, L., Bordi, J. J., Ekelund, J. E., Roth, observatory LISA [60]. As a concluding remark, we note D. C. Precise cassini navigation during solar conjunctions that GRAPE is currently tailored to process cruise-phase through multifrequency plasma calibrations. Journal of data only. Future modifications will extend the platform Guidance, Control, and Dynamics, 2004, 27(2): 251–257. to planetocentric missions, which will benefit from the [7] Reasenberg, R. D., Shapiro, I. I., MacNeil, P. E., extensive treatment of relativistic reference frames and Goldstein, R. B., Breidenthal, J. C., Brenkle, J. P., Cain, their associated spacetime transformations, time scales, D. L., Kaufman, T. M., Komarek, T. A., Zygielbaum, and relativistic multipole structures readily available in A. I. Viking relativity experiment - Verification of signal the vast and extensive literature (see Refs. [8, 61] for retardation by solar gravity. The Astrophysical Journal further details). Letters, 1979, 234: L219. [8] Soffel, M., Klioner, S. A., Petit, G., Wolf, P., Kopeikin, Note: Some findings of this study were presented at the S. M., Bretagnon, P., Brumberg, V. A., Capitaine, 28th International Symposium on Space Flight Dynamics N., Damour, T., Fukushima, T., et al. The IAU (ISSFD-2022-036) [62]. 2000 resolutions for astrometry, celestial mechanics and metrology in the relativistic framework: Explanatory Acknowledgements supplement. The Astronomical Journal, 2003, 126: 2687– The authors acknowledge a helpful discussion with Prof. Sergei Kopeikin regarding the Orbit Determination [9] Lucchesi, D. M., Iafolla, V. The Non-Gravitational Program. JO’L acknowledges support from the Australian Perturbations impact on the BepiColombo Radio Science 12 J. O’Leary, J.P. Barriot Experiment and the key role of the ISA accelerometer: γ in the Cassini experiment. Physics Letters A, 2007, Direct solar radiation and albedo effects. Celestial 367(4–5): 276–280. Mechanics and Dynamical Astronomy, 2006, 96(2): 99– [23] Bertotti, B., Iess, L., Tortora, P. A test of general relativity using radio links with the Cassini spacecraft. Nature, 2003, 425(6956): 374–376. [10] Cappuccio, P., Di Ruscio, A., Iess, L., Mariani, M. [24] Imperi, L., Iess, L. Testing general relativity during the J. BepiColombo gravity and rotation experiment in a cruise phase of the BepiColombo mission to Mercury. In: pseudo drag-free system. In: Proceedings of the AIAA Proceedings of the 2015 IEEE Metrology for Aerospace, Scitech 2020 Forum, 2020: AIAA 2020-1095. Benevento, Italy, 2015: 135–140. [11] Santoli, F., Fiorenza, E., Lefevre, C., Lucchesi, D. [25] Di Stefano, I., Cappuccio, P., Iess, L. The BepiColombo M., Lucente, M., Magnafico, C., Morbidini, A., Peron, solar conjunction experiments revisited. Classical and R., Iafolla, V. ISA, a high sensitivity accelerometer in Quantum Gravity, 2020, 38(5): 055002. the interplanetary space. Space Science Reviews, 2020, [26] Iafolla, V., Fiorenza, E., Lefevre, C., Morbidini, A., 216(8): 145. Nozzoli, S., Peron, R., Persichini, M., Reale, A., Santoli, [12] Iess, L., Asmar, S. W., Cappuccio, P., Cascioli, G., De F. Italian Spring Accelerometer (ISA): A fundamental Marchi, F., di Stefano, I., Genova, A., Ashby, N., Barriot, support to BepiColombo Radio Science Experiments. J. P., Bender, P., et al. Gravity, geodesy and fundamental Planetary and Space Science, 2010, 58(1–2): 300–308. physics with BepiColombo’s MORE investigation. Space [27] Cappuccio, P., Notaro, V., di Ruscio, A., Iess, L., Science Reviews, 2021, 217(1): 21. Genova, A., Durante, D., di Stefano, I., Asmar, S. W., [13] Moyer, T. D. Mathematical formulation of the Double- Ciarcia, S., Simone, L. Report on first inflight data of Precision Orbit Determination Program (DPODP). BepiColombo’s mercury orbiter radio science experiment. Technical Report 32-1527. Jet Propulsion Lab., California IEEE Transactions on Aerospace and Electronic Systems, Inst. Technology, Pasadena, California, USA, 1971. 2020, 56(6): 4984–4988. [14] Moyer, T. D. Formulation for Observed and Computed [28] Cappuccio, P., di Stefano, I., Cascioli, G., Iess, L. Values of Deep Space Network Data Types for Navigation. Comparison of light-time formulations in the post- Hoboken, NJ, USA: John Wiley & Sons, Inc., 2005. Newtonian framework for the BepiColombo MORE [15] Marty, J. Algorithmic documentation of the GINS experiment. Classical and Quantum Gravity, 2021, software. GINS Algorithm Overview, 2013. Available at 38(22): 227001. https://www5.obs-mip.fr/wp-content-omp/uploads/sit [29] Amaro Seoane, P., Arca Sedda, M., Babak, S., Berry, es/28/2017/11/GINS Algo 2013.pdf. C. P. L., Berti, E., Bertone, G., Blas, D., Bogdanovic,´ [16] Evans, S., Taber, W., Drain, T., Smith, J., Wu, H. C., T., Bonetti, M., Breivik, K., et al. The effect of mission Guevara, M., Sunseri, R., Evans, J. MONTE: The next duration on LISA science objectives. General Relativity generation of mission design and navigation software. and Gravitation, 2022, 54(1): 3. CEAS Space Journal, 2018, 10(1): 79–86. [30] Pireaux, S. Time scales in LISA. Classical and Quantum [17] Einstein, A., Infeld, L., Hoffmann, B. The gravitational Gravity, 2007, 24(9): 2271–2281. equations and the problem of motion. The Annals of [31] O’Leary, J., Barriot, J. P. An application of symplectic Mathematics, 1938, 39(1): 65–100. integration for general relativistic planetary orbitography [18] Soffel, M., Langhans, R. Space-Time Reference Systems. subject to non-gravitational forces. Celestial Mechanics Springer Berlin Heidelberg, 2012. and Dynamical Astronomy, 2021, 133(11): 56. [19] Damour, T., Soffel, M., Xu, C. M. General-relativistic [32] Fox, N. J., Velli, M. C., Bale, S. D., Decker, R., Driesman, celestial mechanics. I. Method and definition of reference A., Howard, R. A., Kasper, J. C., Kinnison, J., Kusterer, systems. Physical Review D, 1991, 43(10): 3273–3307. M., Lario, D., et al. The solar probe plus mission: [20] Bertotti, B., Ashby, N., Iess, L. The effect of the motion Humanity’s first visit to our star. Space Science Reviews, of the Sun on the light-time in interplanetary relativity 2016, 204(1): 7–48. experiments. Classical and Quantum Gravity, 2008, [33] Benkhoff, J., Murakami, G., Baumjohann, W., Besse, S., 25(4): 045013. Bunce, E., Casale, M., Cremosese, G., -H Glassmeier, K., [21] Kopeikin, S. M., Sch¨afer, G., Polnarev, A. G., Vlasov, Hayakawa, H., Heyner, D., et al. BepiColombo - Mission I. Y. The orbital motion of Sun and a test of general overview and science goals. Space Science Reviews, 2021, relativity using radio links with the Cassini spacecraft. 217(8): 90. arXiv preprint, 2006, arXiv:gr-qc/0604060. [34] Daquin, J., Alessi, E. M., O’Leary, J., Lemaitre, A., [22] Kopeikin, S. M., Polnarev, A. G., Sch¨afer, G., Vlasov, I. Buzzoni, A. Dynamical properties of the Molniya satellite Y. Gravimagnetic effect of the barycentric motion of the constellation: Long-term evolution of the semi-major Sun and determination of the post-Newtonian parameter axis. Nonlinear Dynamics, 2021, 105(3): 2081–2103. Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 13 [35] Brumberg, V. A. Essential Relativistic Celestial W. M., Park, R. S., Ellis, J. A., Ely, T., Taylor, S. R., Mechanics. CRC Press, 2017. Vallisneri, M. Solar system ephemerides, pulsar timing, [36] Acton, C. H. NASA’s SPICE system models the solar gravitational waves, & navigation. Proceedings of the system. In: Dynamics and Astrometry of Natural and International Astronomical Union, 2017, 13(S337): 150– Artificial Celestial Bodies . Wytrzyszczak, I. M., Lieske, 153. J. H., Feldman, R. A., Eds. Dordrecht: Springer, 1997: [51] Roth, D. C., Guman, M. D., Ionasescu, R. Cassini orbit 257–262. reconstruction from Earth to Jupiter. Root, V1, 2002. [37] Vallado, D. A. Fundamentals of Astrodynamics and Available at https://hdl.handle.net/2014/13042. Applications, Vol. 12. Springer Dordrecht, 2001. [52] Abbate, S. F., Armstrong, J. W., Asmar, S. W., Barbinis, [38] Hill, J. M., O’Leary, J. Generalized transformations E., Bertotti, B., Fleischman, D. U., Gatti, M. S., Goltz, and coordinates for static spherically symmetric general G. L., Herrera, R. G., Iess, L., et al. The Cassini relativity. Royal Society Open Science, 2018, 5(4): gravitational wave experiment. Astronomical Telescopes and Instrumentation. In: Proceedings of the SPIE 4856, [39] Stephani, H., Kramer, D., MacCallum, M., Hoenselaers, Gravitational-Wave Detection, Waikoloa, Hawai’i, USA, C., Herlt, E. Exact Solutions of Einstein’s Field 2003, 4856: 90–97 Equations. Cambridge, UK: Cambridge University Press, [53] Hees, A., Lamine, B., Reynaud, S., Jaekel, M. T., Le Poncin-Lafitte, C., Lainey, V., Fuzfa, ¨ A., Courty, J. M., [40] Pireaux, S., Barriot, J. P., Rosenblatt, P. (SC) RMI: A Dehant, V., Wolf, P. Radioscience simulations in general (S)emi-(C)lassical (R)elativistic (M)otion (I)ntegrator, relativity and in alternative theories of gravity. Classical to model the orbits of space probes around the Earth and and Quantum Gravity, 2012, 29(23): 235027. other planets. Acta Astronautica, 2006, 59(7): 517–523. [54] De Sabbata, V, Melnikov, V. N. Gravitational [41] Burcev, P. Non-gravitational force effect in general theory Measurements, Fundamental Metrology and Constants, of relativity. Cechoslovackij Fiziceskij Zurnal B, 1962, Vol. 230. Springer Science & Business Media, 2012. 12(10): 727–733. [55] Novello, M., Bittencourt, E. Gordon metric revisited. [42] O’Leary, J., Hill, J. M., Bennett, J. C. On the energy Physical Review D, 2012, 86(12): 124024. integral for first post-Newtonian approximation. Celestial [56] Kumar, V. A. A Python-based tool for constructing Mechanics and Dynamical Astronomy, 2018, 130(7): 44. [43] Di Benedetto, M., Iess, L., Roth, D. C. The non- observables from the DSN’s closed-loop archival tracking gravitational accelerations of the Cassini spacecraft. In: data files. SoftwareX, 2022, 19: 101190. Proceedings of the 21st International Symposium on [57] di Stefano, I., Cappuccio, P., Di Benedetto, M., Iess, L. Space Flight Dynamics, Tolouse, 2009. A test of general relativity with ESA’s JUICE mission. [44] Weinberg, S., Wagoner, R. Gravitation and Cosmology: Advances in Space Research, 2022, 70(3): 854–862. Principles and Applications of the General Theory of [58] De Marchi, F., Cascioli, G., Ely, T., Iess, L., Burt, E. Relativity. John Wiley & Sons, 1973. A., Hensley, S., Mazarico, E. Testing the gravitational [45] Hairer, E., Lubich, C., Wanner, G. Geometric redshift with an inner Solar System probe: The VERITAS Numerical Integration: Structure-Preserving Algorithms case. arXiv preprint, 2022, arXiv:2211.08964. for Ordinary Differential Equations. Springer Berlin [59] De Marchi, F., Cascioli, G., Ely, T., Iess, L., Burt, E. Heidelberg, 2006. A., Hensley, S., Mazarico, E. Testing the gravitational [46] Butcher, J. C. Implicit Runge-Kutta processes. redshift with an inner Solar System probe: The VERITAS Mathematics of Computation, 1964, 18(85): 50–64. case. Physical Review D, 2023, 107(6): 064032. [47] O’Leary, J. General Relativistic and Post-Newtonian [60] Danzmann, K., Team, T. L. S. LISA: Laser interferometer Dynamics for Near-Earth Objects and Solar System space antenna for gravitational wave measurements. Bodies. Springer International Publishing, 2021. Classical and Quantum Gravity, 1996, 13(11A): A247. [48] Folkner, W., Williams, J., Boggs, D., Park, R., [61] Klioner, S. A. Relativity in fundamental astronomy: Kuchynka, P. The planetary and lunar ephemerides Solved and unsolved problems. In: Proceedings of the DE430 and DE431. IPN Progress Report 42-196, 2014. Journees ´ Systemes ` de Reference ´ ´ Spatio-temporels, 2007: Available at https://ipnpr.jpl.nasa.gov/progress report/ 42-196/196C.pdf. [49] Fienga, A., Bigot, L., Mary, D., Deram, P., Di Ruscio, [62] O’Leary, J., Barriot, J. A note on the perturbed motion of interplanetary probes in the framework of general A., Bernus, L., Gastineau, M., Laskar, J. Evolution of INPOP planetary ephemerides and bepi-colombo relativity. In: Proceedings of the 28th International simulations. Proc IAU, 2022, 15(S364): 31–51. Symposium on Space Flight Dynamics, Beijing, China, [50] Lazio T Joseph, W., Bhaskaran, S., Cutler, C., Folkner, 2022. 14 J. O’Leary, J.P. Barriot Joseph O’Leary is a postdoctoral interests across a wide range of topics, including geodesy, research fellow at The University of planetary orbitography, geodynamics, and solar physics. E- Melbourne. Dr. O’Leary received his mail: jean-pierre.barriot@upf.pf. Ph.D. degree from the University of South Australia in 2019 and subsequently worked as an industry research fellow with the Open Access This article is licensed under a Creative Space Environment Research Centre and Commons Attribution 4.0 International License, which EOS Space Systems between 2019 and permits use, sharing, adaptation, distribution and 2022. His research interests include nonlinear state estimation, reproduction in any medium or format, as long as you give interplanetary navigation, and general relativity. E-mail: appropriate credit to the original author(s) and the source, joe.oleary@unimelb.edu.au. provide a link to the Creative Commons licence, and indicate if changes were made. Jean-Pierre Barriot is the director The images or other third party material in this article are of the Geodesy Observatory of Tahiti. included in the article’s Creative Commons licence, unless Prof. Barriot received his Ph.D. degree indicated otherwise in a credit line to the material. If material in physics from Montpellier University is not included in the article’s Creative Commons licence and in 1987 and subsequently worked as a your intended use is not permitted by statutory regulation or research engineer with the French Space exceeds the permitted use, you will need to obtain permission Agency (CNES) from 1989 to 2006. He directly from the copyright holder. is currently a Distinguished Professor of To view a copy of this licence, visit http://creativecomm Geophysics at the University of French Polynesia with research ons.org/licenses/by/4.0/. http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Astrodynamics Springer Journals

Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework: An application to the Cassini gravitational wave experiment

Astrodynamics , Volume 7 (3) – Sep 1, 2023

Loading next page...
 
/lp/springer-journals/reconstructing-the-cruise-phase-trajectory-of-deep-space-probes-in-a-olpt03FCK9

References (65)

Publisher
Springer Journals
Copyright
Copyright © The Author(s) 2023
ISSN
2522-008X
eISSN
2522-0098
DOI
10.1007/s42064-023-0160-x
Publisher site
See Article on Publisher Site

Abstract

Astrodynamics https://doi.org/10.1007/s42064-023-0160-x Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework: An application to the Cassini gravitational wave experiment 1,2 3 Joseph O’Leary (B), Jean-Pierre Barriot 1. School of Physics, University of Melbourne, Parkville, VIC 3010, Australia 2. Australian Research Council Centre of Excellence for Gravitational Wave Discovery (OzGrav), Parkville, VIC 3010, Australia 3. Observatoire Geo´ esique d´ de Tahiti, University of French Polynesia, BP 6570, 98702, Faa’a, Tahiti, French Polynesia ABSTRACT KEYWORDS Einstein’s theory of general relativity is playing an increasingly important role in fields general relativity such as interplanetary navigation, astrometry, and metrology. Modern spacecraft and planetary orbitography interplanetary probe prediction and estimation platforms employ a perturbed Newtonian orbit prediction and framework, supplemented with the Einstein–Infeld–Hoffmann n-body equations of motion. determination While time in Newtonian mechanics is formally universal, the accuracy of modern radiometric tracking systems necessitate linear corrections via increasingly complex and error-prone post-Newtonian techniques—to account for light deflection due to the solar system bodies. With flagship projects such as the ESA/JAXA BepiColombo mission now operating at unprecedented levels of accuracy, we believe the standard corrected Newtonian paradigm is approaching its limits in terms of complexity. In this paper, we employ a novel prototype software, General Relativistic Accelerometer-based Propagation Environment, to reconstruct the Cassini cruise-phase trajectory during its first gravitational Research Article wave experiment in a fully relativistic framework. The results presented herein agree with Received: 30 November 2022 post-processed trajectory information obtained from NASA’s SPICE kernels at the order Accepted: 27 February 2023 of centimetres. © The Author(s) 2023 1 Introduction tests of general relativity and facilitated the accurate determination of post-Newtonian parameters [3]. Of The accuracy of interplanetary navigation, planetary particular interest for weak-field tests of gravitation is the sciences, and solar system tests of general relativity has determination of the parameter γ (equal to one in general greatly benefited from modern improvements in multi-link relativity), which measures the amount of space curvature telecommunication configurations and the use of high- produced by a massive body [3]. The deflection and time frequency microwave signals between ground stations delay of photons owing to spacetime curvature can be and interplanetary probes. Precise radiometric tracking expressed in terms of γ , which is currently constrained data obtained during planetary flybys and dedicated −5 to within (2.1 ± 2.3) × 10 of the predicted unity body-centric missions have provided ample opportunity value [4]. By exploiting the multi-frequency transmission for scientific exploration of the mass, gravitational field, available at both ground station and spacecraft, the atmosphere, and orbits of celestial bodies; see Ref. [1] solar conjunction experiment performed by the Cassini and references therein for an overview of notable space probe in 2002 [5] provided the first complete solar plasma missions dedicated to planetary sciences. calibration for both uplink and donwlink signals [6]. The The use of joint Ka/X band microwave links coupled with next-generation multi-frequency technologies [2] Cassini radioscience experiment improved the estimation has yielded unprecedented accuracy for solar system of γ by approximately two orders of magnitude, which B joe.oleary@unimelb.edu.au 2 J. O’Leary, J.P. Barriot −3 was previously constrained to within ∼ 10 of the general [23] correctly account for the relative motion of the relativistic prediction using data from the Viking Martian Sun with respect to the solar system’s barycenter. The lander relativity experiment [7]. expected theoretical difference affects the obtained value −5 Through a series of resolutions, the International of the post-Newtonian parameter γ = (2.1 ± 2.3) × 10 −4 Astronomical Union (IAU) [8] suggests that all problems by the amount ± 1.2 × 10 owing to the barycentric in the field of astronomy or astrodynamics should be velocity of the Sun [2, 21]. Nevertheless, the legacy formulated within the framework of Einstein’s theory software used in the Cassini radioscience experiment, of general relativity. These include: (i) the derivation i.e., Orbit Determination Program [13, 14], employed an of dynamical equations of motion; (ii) the definition approximate model for the motion of the Sun, with a of observables; and (iii) the transformations between negligible effect on the final result [20]. reference frames and time scales. The high precision It is expected that the current best estimate range and Doppler observables combined with data from of γ will be improved in the near future by the onboard accelerometers [9–12] necessitate that modern Mercury Orbiter Radioscience Experiment (MORE) orbitography software employs high-fidelity dynamical onboard the ESA/JAXA BepiColombo mission [24]. and measurement models for spacecraft propagation and The radioscience instrumentation onboard BepiColombo nonlinear state and parameter estimation, respectively will provide plasma-free range information in addition [13, 14]. to stable Ka/X band microwave links [25], while Popular orbitography software such as the French simultaneously benefiting from the direct measurement Space Agency Geo´ desie ´ par Inegrations t´ Numeriques ´ of non-gravitational forces in the local frame of the Simultanees ´ (GINS) [ 15] or NASA’s Mission Analysis, spacecraft by an onboard accelerometer [26]. The Operations, and Navigation Toolkit Environment novel 24 Mcps pseudo-noise ranging instrumentation (MONTE) [16] currently describe the motion of onboard BepiColombo [27] will provide two-way range interplanetary spacecraft via a classical Newtonian information at the centimeter level, so that second- framework linearly corrected with effective forces, order post-Newtonian contributions in the light-time accounting for the effects of general relativity via the so- solution need to be considered in the data reduction called n-body Einstein–Infeld–Hoffmann (EIH) equations process [28]. In addition, future space-based gravitational of motion [17]. Given the stringent accuracy requirements wave detectors, e.g., the Laser Interferometer Space associated with fields such as astrometry, geodesy, Antenna (LISA), are expected to detect the propagation −5 celestial mechanics, and metrological sciences [18], of low-frequency gravitational waves (10 –1 Hz), spacecraft and interplanetary probe orbit propagation complementing the scientific operations of ground- tools based on the so-called “Newton + correction” based detectors such as the Laser Interferometer framework need to include subtle relativistic effects Gravitational wave Observatory (LIGO), VIRGO, and in both state and measurement models to reflect KAGRA collaborations, which focus on higher frequencies the rapid improvements in modern measurement in the gravitational wave frequency spectrum. LISA technologies. The aforementioned ad-hoc process of employs time-delay interferometry (TDI) to reduce reducing relativistic phenomena to simple perturbations the effects of laser frequency and optical bench noise or relativistic corrections introduces a neo-Newtonian [29]. TDI analysis employs solar system barycentric [19] interpretation of the theory of general relativity, coordinate time while each LISA spacecraft records data which can lead to spurious results when equivalence according to the individual spacecraft proper time [30]. between Newtonian and relativistic concepts is assumed. Hence, sophisticated four-dimensional general relativistic Moreover, the “Newton + correction” approach creates spacetime transformations are required for TDI data the likelihood of neglecting potentially important reduction to detect gravitational waves using space-based relativistic contributions if care is not taken. Indeed, it interferometry. Bearing this in mind, we believe that has been discussed at length in relevant literature [20–22] the classical “Newton + correction” paradigm is now whether the dynamical equations and the corresponding reaching its limits in terms of complexity. measurement models used in the data reduction process In a recent article [31], we presented the first results of of the celebrated Cassini solar-conjunction experiments a prototype software, namely General Relativistic Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 3 Accelerometer-based Propagation Environment 2 GRAPE: A brief overview (GRAPE), which describes the perturbed motion In this section, we discuss the mathematical framework of Parker Solar Probe [32], Mercury Planetary Orbiter implemented in GRAPE. First, we provide a high-level [33], and Molniya-like [34] interplanetary probes and overview of the platform and discuss some aspects of spacecraft subject to a radiation-like non-gravitational the software that are unique to general relativity and four-force using the complete framework of general not common with classical propagation environments. relativity. GRAPE employs extended relativistic In Section 2.1, we present the equations of motion and equations of motion that account for non-gravitational discuss the treatment of non-gravitational perturbations forces using either end user supplied accelerometer data in a fully relativistic framework. or approximate dynamical models in the local frame In the framework of general relativity, the unperturbed of the spacecraft. We exploit the tetrad formalism [35] motion of test particles is described using a geodesic to compare dynamical quantities and non-gravitational equation parameterized by an arbitrary affine parameter forces in both co-moving and global frames. The s. For massive test particles such as interplanetary probes novel approach adopted by GRAPE assumes that and near-Earth spacecraft, it is common to parameterize the local structure of spacetime is unperturbed by equations of motion using their associated proper time τ external non-gravitational forces, and that all associated such that s ≡ τ . While proper time parameterization is gravitational and relativistic contributions are contained convenient when working in a theoretical framework, it is and accounted for at the order of the spacetime metric inconvenient for practical purposes. Spacecraft operators, tensor. In this paper, we discuss the preliminary steps tracking stations, and orbit determination as well as data taken to extend GRAPE to the operational planning reduction software reference the state, error covariance, of future missions within the solar system. Further, and observables with respect to a derived external we describe the initial approach adopted to interface coordinate time scale t such as terrestrial, barycentric, it with NASA’s SPICE kernels [36], and verify our or ephemeris times. Hence, GRAPE internally adopts results by comparing the pre-fit Doppler residuals 0 a global coordinate time t = x /c during numerical using GRAPE and SPICE with high-precision Ka/Ka propagation, and the corresponding equations of motion band Doppler data obtained from the first Cassini are transformed from an affine to a coordinate time probe gravitational wave experiment (GWE1) [5]. parameterization accordingly. The remainder of this paper is organized as follows. The incorporation of non-gravitational forces in In Sections 2 and 3, we summarize the equations of classical Newtonian propagation environments is a motion and associated numerical integration procedures well-recognized problem and remains an active area implemented in GRAPE, respectively. In Section 4, we of research in operational astrodynamics [37]. As present new results that extend GRAPE for practical discussed at length in Ref. [31], the treatment of mission planning using radiometric observables and non-gravitational perturbations in general relativity determine that GRAPE recovers the Cassini trajectory has not received the same attention in the literature at the order of centimetres. We conclude the paper in and requires further elucidation. The nonlinear field Section 5, providing a brief discussion on future plans equations of general relativity [38, 39] can be decomposed for the platform. according to contributions from the Einstein tensor Unless otherwise stated, we assume the Einstein G and the energy momentum tensor T . The µν µν summation convention for repeated indices; the use Einstein tensor is composed of the metric tensor of Latin indices (i, j = 1, 2, 3) is reserved for spatial and quantities derived from it, e.g., the Ricci and coordinates, and Greek indices (α, β = 0, 1, 2, 3) denote Riemann tensors, which encode the geometry of spacetime spacetime coordinates, with the 0th component reserved so that all gravitational phenomena are captured α i for time with x = (ct, x ). Finally, we assume a entirely at the order of the spacetime metric tensor. spacetime metric signature according to (+,−, −, − ), Similarly, quantities derived from the energy-momentum a˙ denotes the derivative of a with respect to time, and c tensor encode information related to non-gravitational denotes the speed of light in vacuum. phenomena. In the parlance of operational astrodynamics 4 J. O’Leary, J.P. Barriot i i and interplanetary navigation, gravitational and non- where we introduced the dimensionless quantity β = v /c i i gravitational perturbations are obtained from the metric in terms of the test particle four-velocity v = dx /dt, and and energy-momentum tensors, respectively. Hence, the the Christoffel symbols of the second kind are denoted perturbed equations of motion employed by GRAPE are by Γ . It is worth noting that the components of the βγ derived from the divergence of the energy-momentum four-force in Eq. (5) are given by Eq. (3), and in the tensor, which is linearly decomposed according to appropriate Newtonian limit, Eq. (5) reduces to Newton’s contributions from the energy-momentum of an isolated second law perturbed by an external force f , which can test particle and the corresponding stresses acting on be most easily observed via the explicit calculation of the the system that produce geodesic and perturbed motion, Christoffel symbols given by Eq. (2) in Ref. [42]. On close inspection, Eq. (3) raises two important respectively [40, 41]. practical questions: (i) How are the components of K 2.1 GRAPE: Equations of motion obtained? (ii) Given that non-gravitational forces are To account for non-gravitational forces in Einstein’s measured in the local frame of the spacecraft, how does theory of general relativity, we decompose the energy- GRAPE relate the quantities between local and global momentum conservation equation according to frames, such as the IAU-recommended [8] geocentric and µν µν µν barycentric celestial reference frames? We address both ∇ T ≡ ∇ (Θ + S ) = 0 (1) µ µ questions below. µν where we introduced the stress tensor S , and we The inner product of two four-vectors is a scalar µν interpret Θ as the energy-momentum of an isolated invariant I. In particular, it can be readily shown that µν test particle. In its simplest form, Θ is given in terms in all reference frames, µ µ of mass density ρ and four velocity u ≡ dx /dτ , with µ ν dx dx µν µ ν I ≡ g = c (6) Θ ≡ ρu u . Equation (1) gives rise to two distinct µν dτ dτ µν µν cases: (i) S = 0 and (ii) S ̸= 0. In the absence Differentiating Eq. (6) with respect to proper time and of stresses, it is well established that Eq. (1) yields a substituting the expression given by Eq. (3), we find that µν geodesic equation of motion. In the case where S = ̸ 0, at every point in spacetime, the four-force is orthogonal we find to the four-velocity with µν ν ∇ Θ = ρK (2) µ σ 2 µ f u ≡ K − (K u ) u /c u = 0 (7) µ µ σ µ µν ν where ∇ S = −ρK is the force density associated For Eq. (7) to be satisfied in the local (L) co-moving with an arbitrary non-gravitational four-force f . µ frame of the spacecraft where u = (c, 0, 0, 0), we find Evaluating the term on the left-hand side of Eq. (2) α i that f = 0, K . Hence, in the local co-moving frame, yields [31]: the non-gravitational perturbations are purely spatial and µ ν µ µν µ ν 2 f ≡ u ∇ u = K g − u u /c (3) ν ν are provided to GRAPE by end users using either known dynamical models or data obtained from an onboard where it is clear that geodesic motion is recovered accelerometer. The three-dimensional (Newtonian) non- when K ≡ 0. As previously discussed, Eq. (3) can gravitational forces acting on interplanetary probes be equivalently expressed with respect to a non-affine used in classical orbitography are well known and coordinate time t following several applications of the documented; see Ref. [43] for explicit examples. However, chain rule: −1 α α in the framework of general relativity, the corresponding dx dx dt models for the impulsive and continuous four-forces dt dτ dτ (4) −2 2 α 2 α 2 α are not readily available, posing a theoretical obstacle d x dt d x d t dx = − 2 2 2 for general relativistic orbitography. By modeling the dt dτ dτ dτ dt non-gravitational four-forces in the local frame of the so that the final equations of motion employed by GRAPE spacecraft first, we circumvent this issue, whereupon the are given by 2 i global four-forces can be obtained using the so-called 1 d x i i k i j k 0 i 0 i k = −Γ − 2Γ β − Γ β β + Γ β + 2Γ β β 00 0k jk 00 0k 2 2 tetrad formalism [35]. c dt To compare dynamical quantities between local and 1 dτ 0 i j k 0 i i + Γ β β β − f β − f (5) jk 2 global frames, GRAPE employs the tetrad formalism c dt Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 5 as follows. At every point in spacetime, we construct a admits first integrals or quadratic invariant quantities, pseudo-Cartesian reference system by introducing e symplectic integration procedures ensure that these (ν ) which satisfy remain constant (or bounded in some cases) throughout ν (α) (β ) the numerical integration. δ = g e e , g = e e η (8) αβ µν µν αβ µ ν (α) (β ) The dynamical systems associated with classical We note that e consists of four orthonormal vectors as (ν ) orbitography software are both time- and velocity- opposed to a single tensor quantity known collectively as dependent (e.g., atmospheric drag perturbations), so a tetrad or vierbein [44]. The vector of interest is indicated that the resulting system of equations are not readily by the subscript in parentheses and the component of integrable. Hence, in a classical framework, structure- that vector is denoted by the superscript. preserving integration schemes do not provide an An important feature of the tetrad formalism is advantage over non-symplectic algorithms. This is not that quantities projected onto the tetrad have physical the case for general relativistic orbitography. Given that meaning. Therefore, Eq. (6) is valid in every reference frame, regardless of (0) (0) µ (i) (i) ν dx = e dx , dx = e dx (9) µ ν time- or velocity-dependent external non-gravitational denote locally measured spatial and temporal intervals four-forces, a unique opportunity to exploit structure- [35]. From Eq. (9), it is clear that the force components preserving integration schemes is available with general are obtained according to relativistic mission planning platforms. GRAPE employs (σ ) (σ ) µ ν ν (σ ) f = e f , f = e f (10) the entire suite of numerical integration procedures µ (σ ) presented in Ref. [46], where the reader is referred Equation (10) provides a direct transformation between to Refs. [31, 46, 47] for the specific Butcher tableau local and global force components when e denotes the (ν ) coefficients required for practical implementation. In components of the co-moving tetrad. If, for example, this section, we review implicit Runge–Kutta algorithms e denoted the components of the so-called natural (ν ) and present a derivation of their structure-preserving tetrad, then a further Lorentz transformation is required properties for quadratic invariants. Finally, we note that to account for the relative motion (see Ref. [31, Section this section is based on Ref. [31, Section 4] and remains 2.4] for further details). largely unchanged from its original presentation. 3 GRAPE: Symplectic integration 3.1 GRAPE: Implicit Runge–Kutta proce- scheme dures In general, classical numerical integration algorithms do not preserve the first integrals (constants of Implicit Runge–Kutta algorithms proceed as follows. motion) associated with dynamical systems; thus, over Given an arbitrary system of N differential equations long integration timescales, energy dissipation occurs, with associated initial conditions: leading to incorrect qualitative dynamical behavior. dx i j i i = f (t, x ), x (t ) = x (i, j = 1, 2, . . . , N) Symplectic integration schemes [45] are advanced dt (11) numerical procedures that preserve the symplectic the solution is propagated x → x with an integration n n+1 structure of Hamiltonian systems such that the so-called step-size h according to flow of the system defines a canonical transformation and belongs to a class of numerical methods known i i i x = x + h b k (12) n+1 n j as structure-preserving geometric integrators. In the j=1 celebrated work of Butcher [46], the author introduces i i i i and discusses the properties of implicit Runge–Kutta k = f t + c h, x + h a k (13) n j jl j n l l=1 schemes, presenting the Butcher tableau coefficients for a where Eq. (13) is defined implicitly. Presently, GRAPE series of s-stage, nth-order algorithms, known collectively determines the individual k using fixed-point iteration, as Gauss collocation methods [45]. Implicit Runge–Kutta schemes are symplectic by definition (see Section 3.2 for where the initial value is assumed to be zero everywhere further details). Hence, given a dynamical system that so that (k ) = 0. Improved strategies will be employed j 6 J. O’Leary, J.P. Barriot in future work to initialize k , which requires further where Eq. (20) is obtained by squaring Eq. (12). investigation [45]. Furthermore, x denotes the ith component of an arbitrary vector at time t = t , and k denotes the mth 3.2 Quadratic invariants approximate slope associated with the ith differential equation. For Eq. (16) to be classified as a first integral, General relativity gives rise to a unique dynamical we must have Q = Q = . . . = Q = C,∀n, where quantity (see Eq. (6)) formally known as a quadratic n+1 n 0 invariant [45]. Thus, the symplectic integration Q ≡ Q(x ) denotes the value of the quadratic function procedures implemented in GRAPE provide an additional (16) at time t and C is an arbitrary constant. Hence, layer of verification during the propagation interval. The using Eq. (13), we can equivalently express Eq. (20) coefficients b and a identified in Ref. [ 46] satisfy the according to i ij constraint: i i j i j M x x = M x x + 2h b M X f ij ij ℓ ij n+1 n+1 n n ℓ b a + b a − b b = 0 (i, j = 1, 2, . . . , s) (14) i ij j ji i j ℓ=1 s s X X which consequently indicates that Gauss collocation 2 i j + h (b b − b a − b a ) M k k (21) ℓ m ℓ ℓm m mℓ ij ℓ m methods [46] are symplectic and completely preserve ℓ=1 m=1 quadratic invariants. The proof is outlined below, and i i i where, for brevity, we have introduced k = f (T, X ), the reader is referred to Ref. [45, Chapter 4] for the s i i i with T = t + c h and X = x + h a k . n j jl n n l=1 l complete details. According to Eqs. (17) and (14), the second and final A first integral I(x ) associated with an arbitrary terms on the right-hand side of Eq. (21) are zero; system of differential equations (11) satisfies the condition therefore, for Runge–Kutta methods satisfying Eq. (14), in Eq. (15): j i i j M x x = M x x . Thus, with respect to the ij ij n+1 n+1 n n dI ∂I α α α i j current problem, we find Q(x ) = I(η (x ), u ), where µν = f (t, x ) = 0 (15) dt ∂x I is defined according to Eq. (6), and the condition given We define a quadratic function of Eq. (11) according to by Eq. (19) corresponds to the orthogonality condition k i j Q(x ) = M x x (16) ij associated with the four-velocity and four acceleration; see Eq. (7). It is important to note that in order to enforce where M is an arbitrary symmetric matrix. Hence, ij strict symplecticity, we must express I in the local co- Q(x ) is the first integral of Eq. (11) if Eq. (17) holds moving tetrad so that at each numerical integration step, [45]: the spacetime is assumed to be flat, which is a formal i j k M x f (t, x ) = 0 (17) ij requirement according to the condition dM /dt = 0. ij where it is clear that we must have dM /dt = 0 which ij To conclude this section, we reiterate the important is demonstrated as Eq. (18): point: In classical orbit problems, symplectic integration i j d dM dx dx ij procedures have been shown to have attractive numerical i j i j j i M x x = x x + M x + M x ij ij ij dt dt dt dt properties for the long-term numerical integration of = 0 (18) Hamiltonian systems [45]. Implicit Runge–Kutta schemes and using the symmetry M = M we find [46] are not restricted to Hamiltonian systems. They ij ji are a family of implicitly defined numerical procedures dx i i j M x = M x f = 0 (19) ij ij for the numerical solution of any system of differential dt if and only if M is a constant. equations, e.g., stiff and non-stiff problems [ 45]. The ij According to Eq. (12), we obtain Eq. (20): Gauss collocation methods introduced above are a particular example of an implicit Runge–Kutta scheme j j i i j i M x x = M x x + h b M x k ij ij ℓ ij n+1 n+1 n n n with the added property that their Butcher tableau ℓ=1 coefficients satisfy the constraint given by Eq. (14). i j + h b M k x Hence, they naturally preserve quadratic invariants, as m ij m n m=1 discussed above. s s X X In many classical problems, the introduction of 2 i j + h b b M k k (20) ℓ m ij ℓ m velocity- and time-dependent perturbations results in ℓ=1 m=1 Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 7 ④ ⑤ a system of differential equations that may not readily Doppler ; and (iii) using the SPICE spkezr program, admit first integrals, such that the structure-preserving we perform a direct comparison of the GRAPE ephemeris properties of symplectic integration schemes can no longer data with the Cassini SPICE kernel and project the be exploited. The equivalent problem in the framework residuals onto the radial, along-, and cross-track axes in of general relativity is different; regardless of the specific the local time-dependent spacecraft coordinate system; non-gravitational four-forces entering Eq. (5), the inner see Section 3.3.3 in Ref. [37] for further details. In the product of two four vectors will always be an invariant. In present analysis, the non-gravitational perturbations particular, Eq. (6) is always satisfied, despite the motion entering Eq. (23) are inferred directly from the Cassini of the spacecraft or probe being subject to external forces. SPICE kernel, and we assume that they are held constant throughout the experiment. Given that this is the first Accordingly, we choose to employ an implicit Runge– application of GRAPE to real radiometric tracking data, Kutta scheme for the present analysis, which ensures we persevere with this assumption and note that it is that Eq. (6) does not diverge from c beyond machine an approximation. We further note that the GRAPE precision, providing an additional accuracy check on estimation module is currently being developed, where modifications to the platform. While the results presented specific time-dependent non-gravitational force models in Section 4 are restricted by the duration of GWE1, i.e., will be applied to parameter estimation problems for O (days), we note that GRAPE can be used for long-term probes during interplanetary cruise. The explicit details integration, e.g., generating solar system ephemerides and steps required to reproduce the present analysis are [O (centuries)] [48, 49], which is of interest to the pulsar discussed below. Finally, we note that by construction, timing community for the detection of gravitational the routines utilized internally by the SPICE toolkit waves using time of arrival residuals [50]. The necessary automatically take into account the Earth ephemeris, modifications to GRAPE are minor and will be explored the location of Earth-based tracking stations, and the in a future study. Earth’s rotation. For approximately 40 days starting November 26, 4 Results 2001, the Cassini probe was tracked during solar In this section, we demonstrate the modifications made opposition in search of gravitational waves in the to GRAPE since Ref. [31] in preparation for future millihertz frequency band [52]. Given that gravitational practical mission planning, which we verify through three wave detection methods depend on the timescale of independent strategies: (i) using GRAPE, we generate the radiation, ground-based detectors based on laser an ephemeris for the Cassini probe for the duration of interferometry cannot separate signals from seismic −5 GWE1 and ensure that Eq. (6) is preserved throughout and other noise sources in the frequency range 10 – −2 the integration period; (ii) using the SPICE mkspk 10 Hz. Consequently, searches for low-frequency program, we construct a SPICE kernel containing the gravitational waves must employ space-based detectors GRAPE-generated Cassini trajectory. We calculate the [5]. The gravitational wave experiments performed light-time solution using the GRAPE and Cassini SPICE by Cassini were, at the time, the most sensitive ②③ kernels at each observation epoch and compare the Doppler tracking campaigns achieved using interplanetary corresponding values with the observed GWE1 two-way probes. Moreover, at a distance of approximately ∼ 8 AU, Cassini became the largest gravitational wave ① See https://naif.jpl.nasa.gov/pub/naif/toolkit docs/FORTRA N/ug/mkspk.html for further details. antenna ever used [5]. The success of the Cassini ② See https://naif.jpl.nasa.gov/pub/naif/pds/data/co-s j e v-spi gravitational wave experiments and solar system tests of ce-6-v1.0/cosp 1000/data/spk/ to download the relevant Cassini general relativity was mediated by two major technical SPICE kernel and https://naif.jpl.nasa.gov/pub/naif/pds/data/ co-s j e v-spice-6-v1.0/cosp 1000/data/spk/041014R SCPSE 0 upgrades. First, Cassini was equipped with a novel 1066 04199.lbl for specific kernel details. multi-link telecommunication configuration with two ③ It is important to note that the trajectory information contained in individual Cassini SPICE kernels are not observables. Rather, they are the result of an orbit determination fit using two-way ④ The dataset ID used here is CO-X-RSS-1-GWE1-V1.0. It can be Doppler and range observables obtained from Cassini’s two low downloaded at NASA’s Planetary Data System. gain antennae and high gain antenna; see Ref. [51] for explicit ⑤ See https://naif.jpl.nasa.gov/pub/naif/toolkit docs/FORTRA details on reconstructing the Cassini trajectory. N/spicelib/spkezr.html for further details. 8 J. O’Leary, J.P. Barriot 2 2 2β − 1 µ s˙ s˙ uplink carriers in the X and Ka bands as well as three k i j − + γ + (1 + γ ) c r c c downlink carriers in the X, Ka1 (coherent with the X- jk k̸=j uplink), and Ka2 (coherent with the Ka-uplink) bands, i j j j 2 (1 + γ ) 3 x x˙ − x x˙ i j − x˙ x˙ − considerably reducing the electromagnetic wave phase 2 2 c 2c r ij scintillation owing to solar plasma [6, 52]. Second, the j j i j radioscience experiments performed by Cassini benefited + x x¨ − x x¨ 2c from a dedicated advanced media calibration system 1 µ i j i + x − x (2 + 2γ ) x˙ located at the Goldstone Deep Space Communications 2 3 c r ij j̸=i Complex (DSS-25). The novel system consisted of water 3 + 4γ µ x¨ vapor radiometers, pressure sensors, and microwave j i j − (1 + 2γ ) x˙ x˙ − x˙ + 2c r ij temperature profilers that calibrated the frequency j̸=i modulation caused by the wet and dry components of 1 dτ 0 i i − f β − f (23) the troposphere [5, 6, 52]. c dt The curved spacetime associated with the celestial i i where we introduced the notation x x = x · x to −5 bodies is described (up to order c in the time denote the inner product between two vector quantities, −3 component, up to order c in the spatial components, and we assume that the components of x are given −4 and up to order c in the off-diagonal components) by in terms of Cartesian coordinates. Equation (23) is the components of the n-body metric tensor (see Ref. [14, defined implicitly in terms of x¨ , such that the formal Section 2]): numerical solution is obtained by iteration. However,   X X µ x˙ for the present demonstration, we assume that terms 2 + 2γ j 2γ µ j j   g = , g = −δ 1 + , 0k kl kl i 2 3 2 c r c r of the form x¨ /c can be replaced by their Newtonian ij ij j̸=i j̸=i counterparts. It is important to note that the final term   X X X µ s˙ 2 µ 2β µ 1 + 2γ appearing in Eq. (23) completely captures the general j j j   g = 1− + − 2 4 4 c r c r c r relativistic contributions of an arbitrary non-gravitational ij ij ij j̸=i j̸=i j̸=i four-force f . X X X 2 (2β − 1) µ µ 1 ∂ r j k ij + − µ The calculation of observable quantities in a complete 4 4 2 c r r c ∂t ij jk j̸=i k̸=j j̸=i general relativistic framework [53] relies on highly (22) accurate light-time solutions and is profoundly different where the gravitational parameter of the jth celestial from classical signal-processing paradigms [54, Chapter body is given by µ ; the kth component of its coordinate 2]. In the present study, light-time is calculated using the k 2 k velocity is given by x˙ ; s˙ = x˙ ; β and γ denote j j j approach presented in Ref. [14, Chapter 8]. We note that post-Newtonian parameters; the relative barycentric light-time solutions obtained by numerically solving the distance between bodies i and j is given by r ; and we ij null geodesic equations where the effects of solar plasma have introduced the Kronecker delta function according are contained within a modified metric tensor [ 55] are to δ . In the present demonstration, we assume that kl currently being explored. body i refers to the Cassini probe and bodies j = Figure 1 displays the temporal evolution of the 1, . . . , n − 1 (j ̸= i) refer to the Sun, the nine planets, quadratic invariant I given by Eq. (16). The quadratic and the Moon, where the initial conditions of all n bodies invariant I depends explicitly on the components of the are obtained from the appropriate SPICE kernels on Nov α numerically integrated four-velocity u and implicitly 26, 2001 at 05:00:00 Barycentric Dynamical Time using α on the numerically integrated spacetime coordinates x the SPICE spkezr program. The extended equations of through the metric tensor components g given by µν motion implemented in GRAPE follow directly from the Eq. (22). At each integration step, GRAPE internally definition of the metric tensor components and Eq. (5), checks whether I deviates from machine precision. and are given by Eq. (4-26) in Ref. [14, Section 4]: This is employed as an additional accuracy check j i 2 i X X d x µ x − x 2 (β + γ ) µ l on modifications made to GRAPE. It is clear from = 1 − 2 2 dt r c r il Fig. 1 that the numerically determined four-velocity ij j̸=i l̸=i Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 9 −16 1.0 ×10 SPICE 2.0 GRAPE 0.5 1.5 1.0 0.0 0.5 −0.5 0.0 −1.0 −0.5 −1.0 −1.5 −1.5 −2.0 −2.0 0 5 10 15 20 25 30 35 40 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 Time (day) Time (day) Fig. 1 Time evolution of the quadratic invariant for the Fig. 2 Doppler residuals using SPICE (black) and GRAPE Cassini probe during GWE1. The orbit is propagated using (cyan). The GRAPE ephemeris is propagated using an an integration time-step h = 300 s with an initial condition integration time-step h = 300 s with an initial condition obtained from NASA’s SPICE kernels on Nov 26, 2001 at obtained from NASA’s SPICE kernels on Nov 26, 2001 05:00:00 Barycentric Dynamical Time. The time axis is given at 05:00:00 Barycentric Dynamical Time. The SPICE in days since the initial epoch. We define δI according to ephemerides are pseudo-observations deduced using a least- 2 2 δI ≡ (I − c )/c , where I is determined at each numerical squares fit from radiometric observables. The time axis is integration step and is defined in Eq. (16). GRAPE preserves given in days since the initial epoch. −16 I in the order of 10 . in orbit determination problems, the goal is to obtain components satisfy Eq. (6) at machine precision, which Gaussian residuals with zero mean. In this regard, the −32 can be further reduced to 10 by selecting 128-bit results presented herein are slightly different. In its (quadruple) precision internally at a computational simplest form, orbit determination software is typically cost. The computational demand for double (64-bit decomposed into two modules: prediction and estimation. precision) and quadruple precision is O (seconds) and The goal of the present study is not estimation via least- O (minutes), respectively. We note that a thorough squares or Kalman filter analysis. Instead, this study aims analysis benchmarking the GRAPE platform is deferred to present the first results of the GRAPE platform during for future study and is outside the scope of the present the prediction phase. As noted above, the estimation paper. module of GRAPE is currently under development and Figure 2 compares the pre-fit Doppler residual using will be the topic of a future paper on general relativistic GRAPE and Cassini SPICE trajectory kernels, as state and parameter estimation. discussed in Steps (4)–(6). The details of the Cassini Figure 3 compares the GRAPE and Cassini ˆ ˆ ˆ kernel and observational data used in the present analysis ephemerides in the time-dependent (R, S, W ) spacecraft are presented in Footnotes ② and ④. For clarity, we frame, as discussed in Step (7). Using GRAPE, the plotted only the residuals for the first pass of Cassini reconstructed Cassini trajectory during cruise agrees −4 during GWE1. However, we note that the difference with the SPICE ephemeris at the level of 10 km. in the residuals remains at the ∼ 0.2 Hz level, with The projected residuals in the cross- and along-track no noticeable secular variation during the first 1 × 10 directions, respectively denoted as ∆ W and ∆ S, are observations. The pre-fit Doppler residuals are slightly negligible. The projected residual along the radial biased owing to uncertainties in the initial state of the direction, denoted as ∆ R, grows secularly over the spacecraft obtained from SPICE. The bias is corrected in duration of the Cassini GWE1. During cruise, the the estimation phase, e.g., by employing a (sequential) dominant non-gravitational force is in the radial direction. batch least-squares or Kalman filter algorithm [ 37]. The Given that we assume the non-gravitational forces remain initial-state bias in Fig. 2 is standard and is comparable to constant during GWE1, the projected residual observed that of other modern orbit prediction and determination along the ∆ R axis is expected. We conjecture that software; see Fig. 1 in Ref. [1] for example. Typically, once specific time-dependent non-gravitational four-force δI Residual (Hz) 10 J. O’Leary, J.P. Barriot at each observation time using finite differences. ∆R 0.0004 (5) Define the non-gravitational (NG) perturbations ∆W acting on the probe according to ∆S 0.0003 a ≃ a − a (24) NG S G The GRAPE (G) and SPICE (S) acceleration vectors 0.0002 are respectively defined by 0.0001 a ≃ a , a ≃ a + a (25) G M S M NG where we have introduced the acceleration owing 0.0000 to the n-body metric tensor components only as a . Under the assumption that a ≃ a during M NG 0 0 5 10 15 20 25 30 35 40 Time (day) the short duration cruise phase (< 40 days), where a denotes a constant vector, we employ Eq. (10) Fig. 3 Comparison of SPICE and GRAPE ephemerides’ residuals. The GRAPE ephemeris is propagated using an at each integration time-step and produce a new integration time-step h = 300 s with an initial condition GRAPE-generated kernel to be used with the SPICE obtained from NASA’s SPICE kernels on Nov 26, 2001 software. The ephemerides associated with the new at 05:00:00 Barycentric Dynamical Time. The SPICE kernel now take into account the non-gravitational ephemerides are pseudo-observations deduced using a least- squares fit from radiometric observables. The time axis is forces acting on the probe under the assumption that given in days since the initial epoch. they remain constant throughout the integration procedure. Again, we reiterate that the adopted models are implemented in GRAPE, the ∆ R residual ⑥ approach is an approximation. The implementation will be significantly reduced . of specific time-dependent models in the local tetrad It is important that the results presented in this is postponed to future work. manuscript are easily reproducible. As a result, and in (6) For Fig. 2, we calculate at each observation time the anticipation of a future public release of the GRAPE two-way light-time solution between the spacecraft software, we outline the specific steps required to and Earth antennae (see Ref. [14, Chapter 8] for reproduce Figs. 1–3: explicit details) and compare the calculated Doppler (1) Generate an ephemeris using GRAPE by numerically frequency with the observed value. solving Eq. (23) subject to the condition that α (7) For Fig. 3, we introduce a time-dependent local f ≡ 0. This step produces an ephemeris for the ˆ ˆ ˆ spacecraft coordinate system (R, S, W ) given in Cassini trajectory during the first gravitational wave terms of the probe position r and velocity vector v experiment under the assumption that the motion as Eq. (26) [37]: of the probe is governed entirely by the n-body solar r r × v system metric tensor components given by Eq. (22). ˆ ˆ ˆ ˆ ˆ R = , W = , S = W × R (26) r |r × v| (2) Using the SPICE program mkspk, create a SPICE We project the residual position vector between kernel from the GRAPE-generated ephemeris GRAPE and SPICE, denoted by δ r ≡ r − r , along G S obtained in Step (1). each of the axes in the local spacecraft frame so that (3) Retrieve the approximate probe coordinates at the residuals plotted in Fig. 3 are given explicitly by each observation time relative to the solar system ˆ ˆ ˆ barycenter using both GRAPE-generated and ∆R ≡ δ r · R, ∆W ≡ δ r · W , ∆S ≡ δ r · S (27) Cassini SPICE kernels, using the SPICE program spkezr. 5 Discussions and conclusions (4) Given the SPICE and GRAPE trajectory We presented the first practical results of a novel general information, numerically calculate the acceleration relativistic orbitography platform, GRAPE. While the ⑥ We performed several experiments with h varying between 15 initial results are in good agreement with NASA’s and 300 s. The results of each experiment were comparable, SPICE software, further work is required. It is of with no noticeable improvement in the solution as h decreases. Accordingly, we adopted h = 300 s. utmost importance that GRAPE can accurately provide Projected residual (km) Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 11 light-time solutions for the calculation of radiometric Research Council Centre of Excellence for Gravitational observables in a fully relativistic framework. Rather Wave Discovery (OzGrav) (Grant No. CE170100004). than adopting the standard approach in the literature, JPB was funded by a DAR grant in planetology from the where light-time is approximated by a series of linearly French Space Agency (CNES). superimposed post-Newtonian corrections that gradually Funding note increase in complexity and are error-prone for higher- order corrections, we are currently implementing light- Open Access funding enabled and organized by CAUL time solutions in GRAPE by numerically solving the null and its Member Institutions. geodesic equation given in Ref. [55]. The latter approach guarantees that all relevant relativistic contributions to Declaration of competing interest the dynamics are appropriately accounted for at the order The authors have no competing interests to declare that of the metric tensor and voids the possibility of neglecting are relevant to the content of this article. potentially important terms [20–22]. The motivation for designing a general relativistic References orbit prediction and estimation tool is not to replace [1] Lari, G., Schettino, G., Serra, D., Tommei, G. Orbit the current popular orbitography platforms based on the determination methods for interplanetary missions: corrected Newtonian approach. Rather, it is to provide Development and use of the Orbit14 software. a fully consistent and highly accurate benchmark for Experimental Astronomy, 2022, 53(1): 159–208. current and future platforms. With radioscience data [2] Bertotti, B., Comoretto, G., Iess, L. Doppler tracking preprocessing software now publicly available (see for of spacecraft with multi-frequency links. Astronomy and example Ref. [56]), archival tracking data from a myriad Astrophysics, 1993, 269(1–2): 608–616. of space missions can be easily reduced, thereby providing [3] Will, C. M. Theory and Experiment in Gravitational Physics. Cambridge University Press, 2018. ample opportunity to revisit historical missions and [4] Bertotti, B., Iess, L., Tortora, P. A test of general verify future modifications to the tool. We expect the relativity using radio links with the Cassini spacecraft. current manuscript to be the first in a series of papers Nature, 2003, 425(6956): 374–376. describing new approaches to general relativistic orbit [5] Kliore, A. J., Anderson, J. D., Armstrong, J. W., Asmar, determination (see for example Ref. [47, Chapter 7]) S. W., Hamilton, C. L., Rappaport, N. J., Wahlquist, H. applicable to missions requiring high precision, e.g., D., Ambrosini, R., Flasar, F. M., French, R. G., et al. BepiColombo [12], JUICE [57], and VERITAS [58, 59] Cassini radio science. In: The Cassini-Huygens Mission. missions, as well as the space-based gravitational wave Dordrecht: Kluwer Academic Publishers, 2005: 1–70. [6] Tortora, P., Iess, L., Bordi, J. J., Ekelund, J. E., Roth, observatory LISA [60]. As a concluding remark, we note D. C. Precise cassini navigation during solar conjunctions that GRAPE is currently tailored to process cruise-phase through multifrequency plasma calibrations. Journal of data only. Future modifications will extend the platform Guidance, Control, and Dynamics, 2004, 27(2): 251–257. to planetocentric missions, which will benefit from the [7] Reasenberg, R. D., Shapiro, I. I., MacNeil, P. E., extensive treatment of relativistic reference frames and Goldstein, R. B., Breidenthal, J. C., Brenkle, J. P., Cain, their associated spacetime transformations, time scales, D. L., Kaufman, T. M., Komarek, T. A., Zygielbaum, and relativistic multipole structures readily available in A. I. Viking relativity experiment - Verification of signal the vast and extensive literature (see Refs. [8, 61] for retardation by solar gravity. The Astrophysical Journal further details). Letters, 1979, 234: L219. [8] Soffel, M., Klioner, S. A., Petit, G., Wolf, P., Kopeikin, Note: Some findings of this study were presented at the S. M., Bretagnon, P., Brumberg, V. A., Capitaine, 28th International Symposium on Space Flight Dynamics N., Damour, T., Fukushima, T., et al. The IAU (ISSFD-2022-036) [62]. 2000 resolutions for astrometry, celestial mechanics and metrology in the relativistic framework: Explanatory Acknowledgements supplement. The Astronomical Journal, 2003, 126: 2687– The authors acknowledge a helpful discussion with Prof. Sergei Kopeikin regarding the Orbit Determination [9] Lucchesi, D. M., Iafolla, V. The Non-Gravitational Program. JO’L acknowledges support from the Australian Perturbations impact on the BepiColombo Radio Science 12 J. O’Leary, J.P. Barriot Experiment and the key role of the ISA accelerometer: γ in the Cassini experiment. Physics Letters A, 2007, Direct solar radiation and albedo effects. Celestial 367(4–5): 276–280. Mechanics and Dynamical Astronomy, 2006, 96(2): 99– [23] Bertotti, B., Iess, L., Tortora, P. A test of general relativity using radio links with the Cassini spacecraft. Nature, 2003, 425(6956): 374–376. [10] Cappuccio, P., Di Ruscio, A., Iess, L., Mariani, M. [24] Imperi, L., Iess, L. Testing general relativity during the J. BepiColombo gravity and rotation experiment in a cruise phase of the BepiColombo mission to Mercury. In: pseudo drag-free system. In: Proceedings of the AIAA Proceedings of the 2015 IEEE Metrology for Aerospace, Scitech 2020 Forum, 2020: AIAA 2020-1095. Benevento, Italy, 2015: 135–140. [11] Santoli, F., Fiorenza, E., Lefevre, C., Lucchesi, D. [25] Di Stefano, I., Cappuccio, P., Iess, L. The BepiColombo M., Lucente, M., Magnafico, C., Morbidini, A., Peron, solar conjunction experiments revisited. Classical and R., Iafolla, V. ISA, a high sensitivity accelerometer in Quantum Gravity, 2020, 38(5): 055002. the interplanetary space. Space Science Reviews, 2020, [26] Iafolla, V., Fiorenza, E., Lefevre, C., Morbidini, A., 216(8): 145. Nozzoli, S., Peron, R., Persichini, M., Reale, A., Santoli, [12] Iess, L., Asmar, S. W., Cappuccio, P., Cascioli, G., De F. Italian Spring Accelerometer (ISA): A fundamental Marchi, F., di Stefano, I., Genova, A., Ashby, N., Barriot, support to BepiColombo Radio Science Experiments. J. P., Bender, P., et al. Gravity, geodesy and fundamental Planetary and Space Science, 2010, 58(1–2): 300–308. physics with BepiColombo’s MORE investigation. Space [27] Cappuccio, P., Notaro, V., di Ruscio, A., Iess, L., Science Reviews, 2021, 217(1): 21. Genova, A., Durante, D., di Stefano, I., Asmar, S. W., [13] Moyer, T. D. Mathematical formulation of the Double- Ciarcia, S., Simone, L. Report on first inflight data of Precision Orbit Determination Program (DPODP). BepiColombo’s mercury orbiter radio science experiment. Technical Report 32-1527. Jet Propulsion Lab., California IEEE Transactions on Aerospace and Electronic Systems, Inst. Technology, Pasadena, California, USA, 1971. 2020, 56(6): 4984–4988. [14] Moyer, T. D. Formulation for Observed and Computed [28] Cappuccio, P., di Stefano, I., Cascioli, G., Iess, L. Values of Deep Space Network Data Types for Navigation. Comparison of light-time formulations in the post- Hoboken, NJ, USA: John Wiley & Sons, Inc., 2005. Newtonian framework for the BepiColombo MORE [15] Marty, J. Algorithmic documentation of the GINS experiment. Classical and Quantum Gravity, 2021, software. GINS Algorithm Overview, 2013. Available at 38(22): 227001. https://www5.obs-mip.fr/wp-content-omp/uploads/sit [29] Amaro Seoane, P., Arca Sedda, M., Babak, S., Berry, es/28/2017/11/GINS Algo 2013.pdf. C. P. L., Berti, E., Bertone, G., Blas, D., Bogdanovic,´ [16] Evans, S., Taber, W., Drain, T., Smith, J., Wu, H. C., T., Bonetti, M., Breivik, K., et al. The effect of mission Guevara, M., Sunseri, R., Evans, J. MONTE: The next duration on LISA science objectives. General Relativity generation of mission design and navigation software. and Gravitation, 2022, 54(1): 3. CEAS Space Journal, 2018, 10(1): 79–86. [30] Pireaux, S. Time scales in LISA. Classical and Quantum [17] Einstein, A., Infeld, L., Hoffmann, B. The gravitational Gravity, 2007, 24(9): 2271–2281. equations and the problem of motion. The Annals of [31] O’Leary, J., Barriot, J. P. An application of symplectic Mathematics, 1938, 39(1): 65–100. integration for general relativistic planetary orbitography [18] Soffel, M., Langhans, R. Space-Time Reference Systems. subject to non-gravitational forces. Celestial Mechanics Springer Berlin Heidelberg, 2012. and Dynamical Astronomy, 2021, 133(11): 56. [19] Damour, T., Soffel, M., Xu, C. M. General-relativistic [32] Fox, N. J., Velli, M. C., Bale, S. D., Decker, R., Driesman, celestial mechanics. I. Method and definition of reference A., Howard, R. A., Kasper, J. C., Kinnison, J., Kusterer, systems. Physical Review D, 1991, 43(10): 3273–3307. M., Lario, D., et al. The solar probe plus mission: [20] Bertotti, B., Ashby, N., Iess, L. The effect of the motion Humanity’s first visit to our star. Space Science Reviews, of the Sun on the light-time in interplanetary relativity 2016, 204(1): 7–48. experiments. Classical and Quantum Gravity, 2008, [33] Benkhoff, J., Murakami, G., Baumjohann, W., Besse, S., 25(4): 045013. Bunce, E., Casale, M., Cremosese, G., -H Glassmeier, K., [21] Kopeikin, S. M., Sch¨afer, G., Polnarev, A. G., Vlasov, Hayakawa, H., Heyner, D., et al. BepiColombo - Mission I. Y. The orbital motion of Sun and a test of general overview and science goals. Space Science Reviews, 2021, relativity using radio links with the Cassini spacecraft. 217(8): 90. arXiv preprint, 2006, arXiv:gr-qc/0604060. [34] Daquin, J., Alessi, E. M., O’Leary, J., Lemaitre, A., [22] Kopeikin, S. M., Polnarev, A. G., Sch¨afer, G., Vlasov, I. Buzzoni, A. Dynamical properties of the Molniya satellite Y. Gravimagnetic effect of the barycentric motion of the constellation: Long-term evolution of the semi-major Sun and determination of the post-Newtonian parameter axis. Nonlinear Dynamics, 2021, 105(3): 2081–2103. Reconstructing the cruise-phase trajectory of deep-space probes in a general relativistic framework 13 [35] Brumberg, V. A. Essential Relativistic Celestial W. M., Park, R. S., Ellis, J. A., Ely, T., Taylor, S. R., Mechanics. CRC Press, 2017. Vallisneri, M. Solar system ephemerides, pulsar timing, [36] Acton, C. H. NASA’s SPICE system models the solar gravitational waves, & navigation. Proceedings of the system. In: Dynamics and Astrometry of Natural and International Astronomical Union, 2017, 13(S337): 150– Artificial Celestial Bodies . Wytrzyszczak, I. M., Lieske, 153. J. H., Feldman, R. A., Eds. Dordrecht: Springer, 1997: [51] Roth, D. C., Guman, M. D., Ionasescu, R. Cassini orbit 257–262. reconstruction from Earth to Jupiter. Root, V1, 2002. [37] Vallado, D. A. Fundamentals of Astrodynamics and Available at https://hdl.handle.net/2014/13042. Applications, Vol. 12. Springer Dordrecht, 2001. [52] Abbate, S. F., Armstrong, J. W., Asmar, S. W., Barbinis, [38] Hill, J. M., O’Leary, J. Generalized transformations E., Bertotti, B., Fleischman, D. U., Gatti, M. S., Goltz, and coordinates for static spherically symmetric general G. L., Herrera, R. G., Iess, L., et al. The Cassini relativity. Royal Society Open Science, 2018, 5(4): gravitational wave experiment. Astronomical Telescopes and Instrumentation. In: Proceedings of the SPIE 4856, [39] Stephani, H., Kramer, D., MacCallum, M., Hoenselaers, Gravitational-Wave Detection, Waikoloa, Hawai’i, USA, C., Herlt, E. Exact Solutions of Einstein’s Field 2003, 4856: 90–97 Equations. Cambridge, UK: Cambridge University Press, [53] Hees, A., Lamine, B., Reynaud, S., Jaekel, M. T., Le Poncin-Lafitte, C., Lainey, V., Fuzfa, ¨ A., Courty, J. M., [40] Pireaux, S., Barriot, J. P., Rosenblatt, P. (SC) RMI: A Dehant, V., Wolf, P. Radioscience simulations in general (S)emi-(C)lassical (R)elativistic (M)otion (I)ntegrator, relativity and in alternative theories of gravity. Classical to model the orbits of space probes around the Earth and and Quantum Gravity, 2012, 29(23): 235027. other planets. Acta Astronautica, 2006, 59(7): 517–523. [54] De Sabbata, V, Melnikov, V. N. Gravitational [41] Burcev, P. Non-gravitational force effect in general theory Measurements, Fundamental Metrology and Constants, of relativity. Cechoslovackij Fiziceskij Zurnal B, 1962, Vol. 230. Springer Science & Business Media, 2012. 12(10): 727–733. [55] Novello, M., Bittencourt, E. Gordon metric revisited. [42] O’Leary, J., Hill, J. M., Bennett, J. C. On the energy Physical Review D, 2012, 86(12): 124024. integral for first post-Newtonian approximation. Celestial [56] Kumar, V. A. A Python-based tool for constructing Mechanics and Dynamical Astronomy, 2018, 130(7): 44. [43] Di Benedetto, M., Iess, L., Roth, D. C. The non- observables from the DSN’s closed-loop archival tracking gravitational accelerations of the Cassini spacecraft. In: data files. SoftwareX, 2022, 19: 101190. Proceedings of the 21st International Symposium on [57] di Stefano, I., Cappuccio, P., Di Benedetto, M., Iess, L. Space Flight Dynamics, Tolouse, 2009. A test of general relativity with ESA’s JUICE mission. [44] Weinberg, S., Wagoner, R. Gravitation and Cosmology: Advances in Space Research, 2022, 70(3): 854–862. Principles and Applications of the General Theory of [58] De Marchi, F., Cascioli, G., Ely, T., Iess, L., Burt, E. Relativity. John Wiley & Sons, 1973. A., Hensley, S., Mazarico, E. Testing the gravitational [45] Hairer, E., Lubich, C., Wanner, G. Geometric redshift with an inner Solar System probe: The VERITAS Numerical Integration: Structure-Preserving Algorithms case. arXiv preprint, 2022, arXiv:2211.08964. for Ordinary Differential Equations. Springer Berlin [59] De Marchi, F., Cascioli, G., Ely, T., Iess, L., Burt, E. Heidelberg, 2006. A., Hensley, S., Mazarico, E. Testing the gravitational [46] Butcher, J. C. Implicit Runge-Kutta processes. redshift with an inner Solar System probe: The VERITAS Mathematics of Computation, 1964, 18(85): 50–64. case. Physical Review D, 2023, 107(6): 064032. [47] O’Leary, J. General Relativistic and Post-Newtonian [60] Danzmann, K., Team, T. L. S. LISA: Laser interferometer Dynamics for Near-Earth Objects and Solar System space antenna for gravitational wave measurements. Bodies. Springer International Publishing, 2021. Classical and Quantum Gravity, 1996, 13(11A): A247. [48] Folkner, W., Williams, J., Boggs, D., Park, R., [61] Klioner, S. A. Relativity in fundamental astronomy: Kuchynka, P. The planetary and lunar ephemerides Solved and unsolved problems. In: Proceedings of the DE430 and DE431. IPN Progress Report 42-196, 2014. Journees ´ Systemes ` de Reference ´ ´ Spatio-temporels, 2007: Available at https://ipnpr.jpl.nasa.gov/progress report/ 42-196/196C.pdf. [49] Fienga, A., Bigot, L., Mary, D., Deram, P., Di Ruscio, [62] O’Leary, J., Barriot, J. A note on the perturbed motion of interplanetary probes in the framework of general A., Bernus, L., Gastineau, M., Laskar, J. Evolution of INPOP planetary ephemerides and bepi-colombo relativity. In: Proceedings of the 28th International simulations. Proc IAU, 2022, 15(S364): 31–51. Symposium on Space Flight Dynamics, Beijing, China, [50] Lazio T Joseph, W., Bhaskaran, S., Cutler, C., Folkner, 2022. 14 J. O’Leary, J.P. Barriot Joseph O’Leary is a postdoctoral interests across a wide range of topics, including geodesy, research fellow at The University of planetary orbitography, geodynamics, and solar physics. E- Melbourne. Dr. O’Leary received his mail: jean-pierre.barriot@upf.pf. Ph.D. degree from the University of South Australia in 2019 and subsequently worked as an industry research fellow with the Open Access This article is licensed under a Creative Space Environment Research Centre and Commons Attribution 4.0 International License, which EOS Space Systems between 2019 and permits use, sharing, adaptation, distribution and 2022. His research interests include nonlinear state estimation, reproduction in any medium or format, as long as you give interplanetary navigation, and general relativity. E-mail: appropriate credit to the original author(s) and the source, joe.oleary@unimelb.edu.au. provide a link to the Creative Commons licence, and indicate if changes were made. Jean-Pierre Barriot is the director The images or other third party material in this article are of the Geodesy Observatory of Tahiti. included in the article’s Creative Commons licence, unless Prof. Barriot received his Ph.D. degree indicated otherwise in a credit line to the material. If material in physics from Montpellier University is not included in the article’s Creative Commons licence and in 1987 and subsequently worked as a your intended use is not permitted by statutory regulation or research engineer with the French Space exceeds the permitted use, you will need to obtain permission Agency (CNES) from 1989 to 2006. He directly from the copyright holder. is currently a Distinguished Professor of To view a copy of this licence, visit http://creativecomm Geophysics at the University of French Polynesia with research ons.org/licenses/by/4.0/.

Journal

AstrodynamicsSpringer Journals

Published: Sep 1, 2023

Keywords: general relativity; planetary orbitography; orbit prediction and determination

There are no references for this article.