distorb

Contents

distorb.c

Subroutines that control the integration of the orbital model.

Author

Russell Deitrick (deitrr)

Date

April 28 2015

Functions

void BodyCopyDistOrb(BODY *dest, BODY *src, int iTideModel, int iNumBodies, int iBody)
void InitializeBodyDistOrb(BODY *body, CONTROL *control, UPDATE *update, int iBody, int iModule)
void InitializeUpdateTmpBodyDistOrb(BODY *body, CONTROL *control, UPDATE *update, int iBody)
void ReadDfCrit(BODY *body, CONTROL *control, FILES *files, OPTIONS *options, SYSTEM *system, int iFile)
void ReadInvPlane(BODY *body, CONTROL *control, FILES *files, OPTIONS *options, SYSTEM *system, int iFile)
void ReadOutputLapl(BODY *body, CONTROL *control, FILES *files, OPTIONS *options, SYSTEM *system, int iFile)
void ReadOutputEigen(BODY *body, CONTROL *control, FILES *files, OPTIONS *options, SYSTEM *system, int iFile)
void ReadOverrideMaxEcc(BODY *body, CONTROL *control, FILES *files, OPTIONS *options, SYSTEM *system, int iFile)
void ReadHaltHillStab(BODY *body, CONTROL *control, FILES *files, OPTIONS *options, SYSTEM *system, int iFile)
void ReadHaltCloseEnc(BODY *body, CONTROL *control, FILES *files, OPTIONS *options, SYSTEM *system, int iFile)
void ReadOrbitModel(BODY *body, CONTROL *control, FILES *files, OPTIONS *options, SYSTEM *system, int iFile)
void InitializeOptionsDistOrb(OPTIONS *options, fnReadOption fnRead[])
void ReadOptionsDistOrb(BODY *body, CONTROL *control, FILES *files, OPTIONS *options, SYSTEM *system, fnReadOption fnRead[], int iBody)
double fndCorrectDomain(double angle)
double fndCalcLongP(double dLongA, double dArgP)
double fndCalcLongA(double dLongP, double dArgP)
void VerifyOrbitModel(CONTROL *control, FILES *files, OPTIONS *options)
void VerifyPericenter(BODY *body, CONTROL *control, OPTIONS *options, char cFile[], int iBody, int iVerbose)
void InitializeHeccDistOrbRD4(BODY *body, UPDATE *update, int iBody, int iPert)
void InitializeKeccDistOrbRD4(BODY *body, UPDATE *update, int iBody, int iPert)
void InitializeHeccDistOrbGR(BODY *body, UPDATE *update, int iBody, int iPert)
void InitializeKeccDistOrbGR(BODY *body, UPDATE *update, int iBody, int iPert)
void InitializePincDistOrbRD4(BODY *body, UPDATE *update, int iBody, int iPert)
void InitializeQincDistOrbRD4(BODY *body, UPDATE *update, int iBody, int iPert)
void VerifyPerturbersDistOrbRD4(BODY *body, int iNumBodies, int iBody)
void InitializeHeccDistOrbLL2(BODY *body, UPDATE *update, int iBody, int iPert)
void InitializeKeccDistOrbLL2(BODY *body, UPDATE *update, int iBody, int iPert)
void InitializePincDistOrbLL2(BODY *body, UPDATE *update, int iBody, int iPert)
void InitializeQincDistOrbLL2(BODY *body, UPDATE *update, int iBody, int iPert)
void VerifyPerturbersDistOrbLL2(BODY *body, int iNumBodies, int iBody)
void VerifyGRCorrLL2(BODY *body, int iNumBodies)
void AssignDistOrbDerivatives(BODY *body, EVOLVE *evolve, UPDATE *update, fnUpdateVariable ***fnUpdate, int iBody)
void NullDistOrbDerivatives(BODY *body, EVOLVE *evolve, UPDATE *update, fnUpdateVariable ***fnUpdate, int iBody)
void VerifyDistOrb(BODY *body, CONTROL *control, FILES *files, OPTIONS *options, OUTPUT *output, SYSTEM *system, UPDATE *update, int iBody, int iModule)
void InitializeUpdateDistOrb(BODY *body, UPDATE *update, int iBody)
void FinalizeUpdateHeccDistOrb(BODY *body, UPDATE *update, int *iEqn, int iVar, int iBody, int iFoo)
void FinalizeUpdateKeccDistOrb(BODY *body, UPDATE *update, int *iEqn, int iVar, int iBody, int iFoo)
void FinalizeUpdatePincDistOrb(BODY *body, UPDATE *update, int *iEqn, int iVar, int iBody, int iFoo)
void FinalizeUpdateQincDistOrb(BODY *body, UPDATE *update, int *iEqn, int iVar, int iBody, int iFoo)
void CountHaltsDistOrb(HALT *halt, int *iNumHalts)
void VerifyHaltDistOrb(BODY *body, CONTROL *control, OPTIONS *options, int iBody, int *iHalt)
int fniHaltHillStab(BODY *body, EVOLVE *evolve, HALT *halt, IO *io, UPDATE *update, fnUpdateVariable ***fnUpdate, int iBody)
int fniHaltCloseEnc(BODY *body, EVOLVE *evolve, HALT *halt, IO *io, UPDATE *update, fnUpdateVariable ***fnUpdate, int iBody)
int fbHaltMaxMutualIncDistorb(BODY *body, EVOLVE *evolve, HALT *halt, IO *io, UPDATE *update, fnUpdateVariable ***fnUpdate, int iBody)

Check the maximum allowed mutual inclination when DistOrb is active.

Parameters:
  • body – A pointer to the current BODY instance

  • evolve – A pointer to the integration EVOLVE instance

  • halt – A pointer to the HALT instance

  • io – A pointer to the IO instance

  • update – A pointer to the UPDATE instance

  • iBody – The current index in the BODY instance, irrelevant in this case because mutual inclination is by definition a multi-body variable

Returns:

TRUE if one mutual incliantion in a system is larger than dHaltMaxMutualInc, FALSE if not

int fbCheckMutualIncDistorb(BODY *body, EVOLVE *evolve, HALT *halt, IO *io, UPDATE *update, int iBody)

Check if mutual inclination exceeds the warning limit

Parameters:
  • body – A pointer to the current BODY instance

  • evolve – A pointer to the integration EVOLVE instance

  • halt – A pointer to the HALT instance

  • io – A pointer to the IO instance

  • update – A pointer to the UPDATE instance

  • iBody – The current index in the BODY instance, irrelevant in this case because mutual inclination is by definition a multi-body variable

Returns:

TRUE if one mutual incliantion in a system is larger than dHaltMaxMutualInc, FALSE if not

void WriteEigen(CONTROL *control, SYSTEM *system)
void WriteBodyDEccDtDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodyDSincDtDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodyDLongPDtDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodyDLongADtDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodyDIncDtDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodySinc(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodyPinc(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodyQinc(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodyDHeccDtDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodyDKeccDtDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodyDPincDtDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void WriteBodyDQincDtDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UNITS *units, UPDATE *update, int iBody, double *dTmp, char cUnit[])
void InitializeOutputDistOrb(OUTPUT *output, fnWriteOutput fnWrite[])
void LogOptionsDistOrb(CONTROL *control, FILE *fp)
void LogDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UPDATE *update, fnWriteOutput fnWrite[], FILE *fp)
void LogBodyDistOrb(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system, UPDATE *update, fnWriteOutput fnWrite[], FILE *fp, int iBody)
void AddModuleDistOrb(CONTROL *control, MODULE *module, int iBody, int iModule)
void PropsAuxDistOrb(BODY *body, EVOLVE *evolve, IO *io, UPDATE *update, int iBody)
void ForceBehaviorDistOrb(BODY *body, MODULE *module, EVOLVE *evolve, IO *io, SYSTEM *system, UPDATE *update, fnUpdateVariable ***fnUpdate, int iBody, int iModule)
unsigned long int fniFactorial(unsigned int n)
int fniNchoosek(int N, int k)

Number of combinations of k in N

Parameters:
  • N – Size of set

  • k – Size of desired subset

Returns:

Binomial coefficient N choose k

int fniCombCount(int x, int y, int N)

Gives the index of a pair of values in N choose 2. For example, for 4 planets, the index for the pair (1,2) -> 0, (1,3) -> 1, (1,4) -> 2, (2,3) -> 3, etc.

Parameters:
  • x – Index of first planet

  • y – Index of second planet

  • N – Number of planets in system

Returns:

Index corresponding to planet pair (x,y)

double fndABmatrix(BODY *body, int j, int jBody, int kBody)

Calculates components AB matrices to find eigenvalues in LL2 problem

Parameters:
  • body – Struct containing all body information and variables

  • j – Index j in Laplace coefficients

  • jBody – Index of perturbed body

  • kBody – Index of perturbing body

Returns:

Component of A or B matrix in rad/year

double fndMutualHillRad(BODY *body, int iBody, int jBody)

Calculates mutual hill radii between two bodies

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

Mutual hill radii in meters

double fndGRCorrMatrix(BODY *body, int jBody, int kBody)

Post-Newtonian correction to AB matrix in LL2 solution

Parameters:
  • body – Struct containing all body information and variables

  • jBody – Index of perturbed body

  • kBody – Index of perturbing body

Returns:

Correction to component of AB matrix (added to A and B) in rad/year

void HessEigen(double **amat, int origsize, double real[], double imag[])

Finds all eigenvalues of an upper Hess. matrix amat

Parameters:
  • amat – Matrix to find eigenvalues of

  • origsize – Size of original matrix

  • real – The real components of the eigenvalues

  • imag – The imaginary components of the eigenvalues (usually ~ 0)

void RowSwap(double **matrix, int size, int i, int j)

Swaps two rows in a matrix

Parameters:
  • matrix – Matrix in question

  • size – The number of rows/columns in matrix (square)

  • i – One of the rows to be swapped

  • j – The other row to be swapped

void ColSwap(double **matrix, int size, int i, int j)

Swaps two columns in a matrix

Parameters:
  • matrix – Matrix in question

  • size – The number of rows/columns in matrix (square)

  • i – One of the columns to be swapped

  • j – The other column to be swapped

void HessReduce(double **a, int size)

Reduces a matrix to Upper Hessenberg form

Parameters:
  • a – Matrix in question

  • size – The number of rows/columns in matrix a (square)

void BalanceMatrix(double **a, int size)

Balances a matrix

Parameters:
  • a – Matrix to be balanced

  • size – The number of rows/columns in matrix a (square)

void LUDecomp(double **amat, double **copy, double *scale, int *rowswap, int size)

Decomposes matrix to LU form

Parameters:
  • amat – Matrix to be LU’ed

  • copy – Copy of matrix containing LU decomposition

  • scale – Vector of scale factors used in decomposition

  • rowswap – Indices of swapped rows

  • size – Size of matrix (square)

void LUSolve(double **lumat, double *soln, int *swap, int size)

Solves system of equations involving an LU matrix

Parameters:
  • lumat – LU matrix

  • soln – Vector containing output solution

  • swap – Indices of swapped rows

  • size – Size of matrix (square)

void FindEigenVecEcc(SYSTEM *system, int count, int pls)

Find eccentricity eigenvectors for LL2 solution

Parameters:
  • system – Struct containing system information

  • count – Index of eigenvalue

  • pls – Number of planets

void FindEigenVecInc(SYSTEM *system, int count, int pls)

Find inclination eigenvectors for LL2 solution

Parameters:
  • system – Struct containing system information

  • count – Index of eigenvalue

  • pls – Number of planets

void SolveEigenVal(BODY *body, EVOLVE *evolve, SYSTEM *system)

Find eigenvalues for LL2 solution

Parameters:
  • body – Struct containing all body information and variables

  • evolve – Struct containing evolve information

  • system – Struct containing system information

void ScaleEigenVec(BODY *body, EVOLVE *evolve, SYSTEM *system)

Scales eigenvectors to initial conditions

Parameters:
  • body – Struct containing all body information and variables

  • evolve – Struct containing evolve information

  • system – Struct containing system information

void RecalcLaplace(BODY *body, EVOLVE *evolve, SYSTEM *system, int iVerbose)

Recalculates Semi-major axis terms in case where RD4 solution is coupled to eqtide

Parameters:
  • body – Struct containing all body information and variables

  • evolve – Struct containing evolve information

  • system – Struct containing system information

  • iVerbose – Verbosity level of output (currently not in use)

void RecalcEigenVals(BODY *body, EVOLVE *evolve, SYSTEM *system)

Recalculates eigenvalues in case where LL2 solution is coupled to eqtide

Parameters:
  • body – Struct containing all body information and variables

  • evolve – Struct containing evolve information

  • system – Struct containing system information

double fndXangle1(BODY *body, int iBody)

First x-term associated rotation into 3D Cartesian coordinates

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of body in question

Returns:

First x-term in rotation into 3D coordinates

double fndXangle2(BODY *body, int iBody)

Second x-term associated rotation into 3D Cartesian coordinates

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of body in question

Returns:

Second x-term in rotation into 3D coordinates

double fndYangle1(BODY *body, int iBody)

First y-term associated rotation into 3D Cartesian coordinates

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of body in question

Returns:

First y-term in rotation into 3D coordinates

double fndYangle2(BODY *body, int iBody)

Second y-term associated rotation into 3D Cartesian coordinates

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of body in question

Returns:

Second y-term in rotation into 3D coordinates

double fndZangle1(BODY *body, int iBody)

First z-term associated rotation into 3D Cartesian coordinates

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of body in question

Returns:

First z-term in rotation into 3D coordinates

double fndZangle2(BODY *body, int iBody)

Second z-term associated rotation into 3D Cartesian coordinates

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of body in question

Returns:

Second z-term in rotation into 3D coordinates

double fndXinit(BODY *body, int iBody)

Calculates x-position in orbital plane

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of body in question

Returns:

Position of planet in x direction (au)

double fndYinit(BODY *body, int iBody)

Calculates y-position in orbital plane

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of body in question

Returns:

Position of planet in y direction (au/day)

double fndVxi(BODY *body, int iBody)

Calculates x-velocity in orbital plane

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of body in question

Returns:

Velocity of planet in x direction (au/day)

double fndVyi(BODY *body, int iBody)

Calculates y-velocity in orbital plane

Parameters:
  • body – Struct containing all body information and variables

  • iBody – Index of body in question

Returns:

Velocity of planet in y direction (au/day)

void cross(double *a, double *b, double *c)

Calculates cross product of vectors

Parameters:
  • a – First vector of cross prodect

  • b – Second vector of cross product

  • c – Resulting product containing cross product

double fndLaplaceCoeff(double dAxRatio, int iIndexJ, double dIndexS)

Laplace coefficient used in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

  • dIndexS – Index s of Laplace Coefficients (usually s = 1/2, 3/2, or 5/2)

Returns:

Laplace coefficient

double fndDerivLaplaceCoeff(int iNthDeriv, double dAxRatio, int iIndexJ, double dIndexS)

Derivative in d/d(alpha) of Laplace coefficient used in disturbing function

Parameters:
  • iNthDeriv – Order of derivative to be taken

  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

  • dIndexS – Index s of Laplace Coefficients (usually s = 1/2, 3/2, or 5/2)

Returns:

Laplace coefficient

double fndSemiMajAxF1(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF1Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF2(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF2Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF3(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF3Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF4(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF4Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF5(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF5Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF6(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF6Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF7(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF7Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF8(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF8Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF9(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF9Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF10(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF10Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF11(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF11Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF12(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF12Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF13(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF13Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF14(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF14Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF15(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF15Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF16(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF16Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF17(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF17Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF18(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF18Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF19(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF19Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF20(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF20Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF21(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF21Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF22(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF22Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF23(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF23Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF24(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF24Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF25(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF25Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndSemiMajAxF26(double dAxRatio, int iIndexJ)

Semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Semi-major axis term

double fndDSemiF26Dalpha(double dAxRatio, int iIndexJ)

Derivative in d/d(alpha) of semi-major axis term in disturbing function

Parameters:
  • dAxRatio – Ratio of inner planet’s semi to outer planet’s (must be < 1)

  • iIndexJ – Index j of Laplace Coefficients (j = 0 for secular model)

Returns:

Derivative of semi-major axis term

double fndDdistDhDir01(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir02(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir03(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir04(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir05(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir06(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir08(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir09(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir010(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir011(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir013(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdistDhDir014(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dh term for interior body

double fndDdisturbDHecc(BODY *body, SYSTEM *system, int *iaBody)

Sums over secular dR/dh terms of disturbing function for interior body

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Total dR/dh for interior body

double fndDdistDhPrmDir01(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir02(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir03(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir04(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir06(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir07(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir09(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir010(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir011(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir012(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir014(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdistDhPrmDir015(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dh of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dh term for exterior body

double fndDdisturbDHeccPrime(BODY *body, SYSTEM *system, int *iaBody)

Sums over secular dR/dh terms of disturbing function for exterior body

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Total dR/dh for exterior body

double fndDdistDkDir01(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir02(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir03(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir04(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir05(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir06(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir08(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir09(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir010(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir011(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir013(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdistDkDir014(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dk term for interior body

double fndDdisturbDKecc(BODY *body, SYSTEM *system, int *iaBody)

Sums over secular dR/dk terms of disturbing function for interior body

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Total dR/dk for interior body

double fndDdistDkPrmDir01(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir02(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir03(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir04(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir06(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir07(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir09(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir010(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir011(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir012(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir014(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdistDkPrmDir015(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dk of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dk term for exterior body

double fndDdisturbDKeccPrime(BODY *body, SYSTEM *system, int *iaBody)

Sums over secular dR/dk terms of disturbing function for exterior body

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Total dR/dk for exterior body

double fndDdistDpDir01(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir02(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir03(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir05(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir06(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir07(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir08(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir09(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir010(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir011(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir012(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdistDpDir016(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dp term for interior body

double fndDdisturbDPinc(BODY *body, SYSTEM *system, int *iaBody)

Sums over secular dR/dp terms of disturbing function for interior body

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Total dR/dp for interior body

double fndDdistDpPrmDir01(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir02(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir03(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir08(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir09(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir010(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir011(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir012(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir013(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir014(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir015(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdistDpPrmDir016(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dp of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dp term for exterior body

double fndDdisturbDPincPrime(BODY *body, SYSTEM *system, int *iaBody)

Sums over secular dR/dp terms of disturbing function for exterior body

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Total dR/dp for exterior body

double fndDdistDqDir01(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir02(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir03(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir05(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir06(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir07(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir08(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir09(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir010(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir011(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir012(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdistDqDir016(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of interior body

  • jBody – Index of exterior body

Returns:

dR/dq term for interior body

double fndDdisturbDQinc(BODY *body, SYSTEM *system, int *iaBody)

Sums over secular dR/dq terms of disturbing function for interior body

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Total dR/dq for interior body

double fndDdistDqPrmDir01(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir02(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir03(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir08(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir09(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir010(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir011(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir012(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir013(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir014(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir015(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdistDqPrmDir016(BODY *body, SYSTEM *system, int iBody, int jBody)

Derivative in d/dq of disturbing function term

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iBody – Index of exterior body

  • jBody – Index of interior body

Returns:

dR/dq term for exterior body

double fndDdisturbDQincPrime(BODY *body, SYSTEM *system, int *iaBody)

Sums over secular dR/dq terms of disturbing function for exterior body

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Total dR/dq for exterior body

double fndApsidalGRCorrection(BODY *body, int *iaBody)

Relativistic precession of periastron

Parameters:
  • body – Struct containing all body information and variables

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative d(longp)/dt due to general relativity

double fndApsidalGRDhDt(BODY *body, SYSTEM *system, int *iaBody)

Relativistic correction to derivative of variable Hecc = e*sin(longp) in RD4 solution

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative dh/dt due to general relativity

double fndApsidalGRDkDt(BODY *body, SYSTEM *system, int *iaBody)

Relativistic correction to derivative of variable Kecc = e*sin(longp) in RD4 solution

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative dk/dt due to general relativity

double fndDistOrbRD4DhDt(BODY *body, SYSTEM *system, int *iaBody)

Derivative of variable Hecc = e*sin(longp) in RD4 solution

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative dh/dt

double fndDistOrbRD4DkDt(BODY *body, SYSTEM *system, int *iaBody)

Derivative of variable Kecc = e*cos(longp) in RD4 solution

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative dk/dt

double fndDistOrbRD4DpDt(BODY *body, SYSTEM *system, int *iaBody)

Derivative of variable Pinc = sin(inc/2)*sin(longa) in RD4 solution

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative dp/dt

double fndDistOrbRD4DqDt(BODY *body, SYSTEM *system, int *iaBody)

Derivative of variable Qinc = sin(inc/2)*cos(longa) in RD4 solution

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative dq/dt

double fndDistOrbLL2Hecc(BODY *body, SYSTEM *system, int *iaBody)

Value of variable Hecc = e*sin(longp) at time dAge, in the LL2 solution (not a derivative)

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Hecc at time dAge

double fndDistOrbLL2Kecc(BODY *body, SYSTEM *system, int *iaBody)

Value of variable Kecc = e*cos(longp) at time dAge, in the LL2 solution (not a derivative)

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Kecc at time dAge

double fndDistOrbLL2Pinc(BODY *body, SYSTEM *system, int *iaBody)

Value of variable Pinc = sin(inc/2)*sin(longp) at time dAge, in the LL2 solution (not a derivative)

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Pinc at time dAge

double fndDistOrbLL2Qinc(BODY *body, SYSTEM *system, int *iaBody)

Value of variable Qinc = sin(inc/2)*cos(longp) at time dAge, in the LL2 solution (not a derivative)

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Qinc at time dAge

double fndDistOrbLL2DhDt(BODY *body, SYSTEM *system, int *iaBody)

Provides derivative of variable Hecc = e*sin(longp) to couple eqtide to LL2 solution.

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative dh/dt

double fndDistOrbLL2DkDt(BODY *body, SYSTEM *system, int *iaBody)

Provides derivative of variable Kecc = e*sin(longp) to couple eqtide to LL2 solution.

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative dk/dt

double fndDistOrbLL2DpDt(BODY *body, SYSTEM *system, int *iaBody)

Provides derivative of variable Pinc = sin(inc/2)*sin(longa) to couple distrot to LL2 solution.

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative dp/dt

double fndDistOrbLL2DqDt(BODY *body, SYSTEM *system, int *iaBody)

Provides derivative of variable Qinc = sin(inc/2)*cos(longa) to couple distrot to LL2 solution.

Parameters:
  • body – Struct containing all body information and variables

  • system – Struct containing system information

  • iaBody – Array containing indices of bodies associated with interaction

Returns:

Derivative dq/dt

double fdInclination(BODY *body, int iBody)
double fdLongA(BODY *body, int iBody)

distorb.h

Subroutines that control the integration of the orbital model.

Author

Russell Deitrick (deitrr)

Date

April 28 2015

Defines

LL2
RD4
TEENY
MAXECCDISTORB
MAXMUTUALINCRD4
MAXMUTUALINCLL2
RADIX
SWAP(g, h)
A(iIndexJ)
B(iIndexJ)
C(iIndexJ)
LAPLNUM
OPTSTARTDISTORB
OPTENDDISTORB
OPT_DFCRIT
OPT_INVPLANE
OPT_ORMAXECC
OPT_HALTHILLSTAB
OPT_HALTCLOSEENC
OPT_EIGENSET
OPT_EIGENVALUE
OPT_EIGENVECTOR
OPT_OUTPUTEIGEN
OPT_OUTPUTLAPL
OPT_ORBITMODEL
OUTSTARTDISTORB
OUTENDDISTORB
OUTBODYSTARTDISTORB
OUT_SINC
OUT_PINC
OUT_QINC
OUT_DECCDTDISTORB
OUT_DINCDTDISTORB
OUT_DSINCDTDISTORB
OUT_DLONGADTDISTORB
OUT_DLONGPDTDISTORB
OUT_DHECCDTDISTORB
OUT_DKECCDTDISTORB
OUT_DPINCDTDISTORB
OUT_DQINCDTDISTORB