Page - 37 - in Contributions to GRACE Gravity Field Recovery - Improvements in Dynamic Orbit Integration, Stochastic Modelling of the Antenna Offset Correction, and Co-Estimation of Satellite Orientations

5.2.3 State Transition Matrix In the first derivations of eq. (5.2.7), the approximate state transition matrix Φ¯ is determined by taking the partial of the integrated velocity and position with regard to the initial state of the orbit arcy0. This approximation neglects that that the terms r˙inte andrinte also depend on the position along the arc, as they are a result of the integration of theaccelerations r¨e evaluatedfromtheforcemodelsat thepositionsre.Thepositions along the arc in turn depend on the initial state, requiring application of the chain rule in the determination ofΦ. Again dividing the problem into position and velocity components, the complete state transition matrix for the positionΦr is arrived at by taking the derivative of eq. (5.1.8) with regard to the initial statey0, giving ∂r(τ) ∂y0 = ∂(r0+τT · r˙0) ∂y0 +T2 ∫ τ 0 (τ−τ′)∂f(τ ′) ∂y0 dτ′ . (5.2.27) Application of the chain rule to the derivative of the force function with regard to the initial state yields ∂f(τ′) ∂y0 = ∂f(τ′) ∂r(τ′) ∂r(τ′) ∂y0 (5.2.28) and thus gives, by inserting into eq. (5.2.27), ∂r(τ) ∂yˆ0︸ ︷︷ ︸ Φr = ∂(r0+τT · r˙0) ∂y0︸ ︷︷ ︸ Φ¯r +T2 ∫ τ 0 (τ−τ′)︸ ︷︷ ︸ Kr ∂f(τ′) ∂r(τ′)︸ ︷︷ ︸ T ∂r(τ′) ∂y0︸ ︷︷ ︸ Φr dτ′ . (5.2.29) Here, the polynomial integration matrixKr, the matrix of Marussi TensorsT , and the initial approximation of the state transition matrix Φ¯r can be identified. In addition, the complete state transition matrixΦ also appears on both sides of the equation system. This system can be solved for the complete state transition matrixΦwith Φr=[I−KrT]−1Φ¯r (5.2.30) where the same inverse as previously encountered in eq. (5.2.25) appears. In fact, both eq. (5.2.25) and eq. (5.2.30) are of the form l=[I−KrT]−1x , (5.2.31) and can thus be solved in a similar fashion. This inverse is of size 3N×3N, containing three position components per epoch. Direct inversion of this matrix is expensive even for moderate arc lengths of a few hours. For arc lengths of 24h at a sampling of 5s, as used in ITSG-Grace2016, this inverse alone would be responsible for a large fraction of the computation time in determining the variational equations. The special blocked structure of the inverse can be exploited to solve the system using efficient algorithms. As integration up to timeτ only depends on accelerations before that point in time, the discretised integration matrixKr is only populated on or 5.2 Orbit Integration and State Transition Matrix 37