interpolate.h#
Functions
-
static inline void comp_sum(real num, real *sum, real *compCoeff)#
Compute the sum of two real numbers with a compensated summation.
- Parameters:
num – [in] Number to add to the sum.
sum – [inout] Sum of the numbers.
compCoeff – [inout] Compensation coefficient.
-
void approx_xInteg_math(const std::vector<real> &xInteg0, const std::vector<real> &accInteg0, const real &dt, const real &h, const std::vector<real> &b, const size_t &dim, const size_t starti, const size_t startb, const size_t &iterStep, std::vector<real> &xIntegNext, std::vector<real> &xIntegCompCoeffs)#
Evaluate the Gauss-Radau polynomial.
- Parameters:
xInteg0 – [in] Initial state vector.
accInteg0 – [in] Initial acceleration vector.
dt – [in] Integration time step.
h – [in] Fraction of the time step to use for the approximation.
b – [in] Interpolation coefficients for the Gauss-Radau polynomial.
dim – [in] Dimension of the system (number of 2nd derivatives).
starti – [in] Starting index of the state vector to approximate.
startb – [in] Starting index of the interpolation coefficients.
iterStep – [in] Number of derivatives to evaluate.
xIntegNext – [out] Approximated state vector.
xIntegCompCoeffs – [out] Compensation coefficients.
-
void approx_xInteg(const std::vector<real> &xInteg0, const std::vector<real> &accInteg0, const real &dt, const real &h, const std::vector<real> &b, const size_t &dim, const std::vector<IntegBody> &integBodies, std::vector<real> &xIntegNext, std::vector<real> &xIntegCompCoeffs)#
Use the Gauss-Radau polynomial to approximate the integral of the acceleration.
- Parameters:
xInteg0 – [in] Initial state vector.
accInteg0 – [in] Initial acceleration vector.
dt – [in] Integration time step.
h – [in] Fraction of the time step to use for the approximation.
b – [in] Interpolation coefficients for the Gauss-Radau polynomial.
dim – [in] Dimension of the system (number of 2nd derivatives).
integBodies – [in] List of integrated bodies in the PropSimulation.
xIntegNext – [out] Approximated state vector.
xIntegCompCoeffs – [out] Compensation coefficients.
-
void interpolate_on_the_fly(PropSimulation *propSim, const real &t, const real &dt)#
Interpolate the integrator state for the time step that was just completed.
- Parameters:
propSim – [in] PropSimulation object for the integration.
t – [in] Time at the beginning of the time step.
dt – [in] Completed time step size.
-
void get_interpIdxInWindow(const PropSimulation *propSim, const real &tWindowStart, const real &tNext, const bool &forwardProp, const bool &backwardProp, bool &interpIdxInWindow)#
Determine whether the next interpolation index is within the window of the time step that was just completed.
- Parameters:
propSim – [in] PropSimulation object for the integration.
tWindowStart – [in] Start time of the interpolation window.
tNext – [in] Time of the next interpolation.
forwardProp – [in] Flag to indicate forward integration.
backwardProp – [in] Flag to indicate backward integration.
interpIdxInWindow – [out] Flag to indicate whether the next interpolation index is within the window.
-
void get_lightTime_and_xRelative(PropSimulation *propSim, const size_t interpIdx, const real tInterpGeom, const std::vector<real> &xInterpGeom, std::vector<real> &lightTime, std::vector<real> &xInterpApparent)#
Compute the light time and apparent position of the target body.
- Parameters:
propSim – [in] PropSimulation object for the integration.
interpIdx – [in] Index of the next interpolation time.
tInterpGeom – [in] Time to interpolate to.
xInterpGeom – [in] Geometric state vector of the target body.
lightTime – [out] Light time to the target body.
xInterpApparent – [out] Apparent state vector of the target body.
-
void get_lightTimeOneBody(PropSimulation *propSim, const size_t &i, const real tInterpGeom, std::vector<real> xInterpGeom, std::vector<real> xObserver, const bool bouncePointAtCenterOfMass, real &lightTimeOneBody)#
Compute the light time to the target body.
- Parameters:
propSim – [in] PropSimulation object for the integration.
i – [in] Index of the target body.
tInterpGeom – [in] Time to interpolate to.
xInterpGeom – [in] Geometric state vector of the target body.
xObserver – [in] State vector of the observer.
bouncePointAtCenterOfMass – [in] Flag to indicate whether the bounce point is at the center of mass (as opposed to leading edge).
lightTimeOneBody – [out] Light time to the target body.