Hi,
I'm using JSBSim as the physics engine for simulating space flight in Outerra Anteworld.
In one of the XML files pertaining to the Jupiter 246 / Ares IV launch vehicle, there are a few equations for calculating the semi-major axis, eccentricity, apoapsis, periapsis and orbital period. I would like to complete the calculations for the rest of the orbital elements, but noticed a serious flaw in the method used for eccentricity, which brings into question the rest of the equations used.
At the start of the simulation, where the rocket is on the launch pad, the eccentricity value is initially 1.0, and during ascent decreases slowly towards 0.0. Immediately I knew something was wrong, as the eccentricity should begin at 0.0 and increase during ascent, until you perform a burn to circularize the orbit, then it should return to near 0.0. I've included the XML equations, as well as what I have re-written into JavaScript functions. They are deriving the Keplerian elements from ECI Cartesian coordinates.
Can anyone determine where they went wrong?
Best regards,
Uriah George
XML:
JavaScript:
I'm using JSBSim as the physics engine for simulating space flight in Outerra Anteworld.
In one of the XML files pertaining to the Jupiter 246 / Ares IV launch vehicle, there are a few equations for calculating the semi-major axis, eccentricity, apoapsis, periapsis and orbital period. I would like to complete the calculations for the rest of the orbital elements, but noticed a serious flaw in the method used for eccentricity, which brings into question the rest of the equations used.
At the start of the simulation, where the rocket is on the launch pad, the eccentricity value is initially 1.0, and during ascent decreases slowly towards 0.0. Immediately I knew something was wrong, as the eccentricity should begin at 0.0 and increase during ascent, until you perform a burn to circularize the orbit, then it should return to near 0.0. I've included the XML equations, as well as what I have re-written into JavaScript functions. They are deriving the Keplerian elements from ECI Cartesian coordinates.
Can anyone determine where they went wrong?
Best regards,
Uriah George
XML:
Code:
<property value="14.07644180e15"> guidance/executive/gm </property>
<channel name="Orbital Parameters">
<fcs_function name="guidance/executive/semimajor-axis-length-ft">
<description>
This equation calculates the semimajor-axis length in feet. a = mu/(2*mu/r - v^2)
</description>
<function>
<quotient>
<p> guidance/executive/gm </p>
<difference>
<quotient>
<product>
<v> 2.0 </v>
<p> guidance/executive/gm </p>
</product>
<p> position/radius-to-vehicle-ft </p>
</quotient>
<pow>
<p> velocities/eci-velocity-mag-fps </p>
<v> 2.0 </v>
</pow>
</difference>
</quotient>
</function>
</fcs_function>
<fcs_function name="guidance/executive/eccentricity">
<description>
This equation calculates eccentricity. e = sqrt (1 - [ r v cos(gamma) ]^2 / {a * mu})
</description>
<function>
<pow>
<difference>
<v> 1.0 </v>
<quotient>
<pow>
<product>
<property> position/radius-to-vehicle-ft </property>
<property> velocities/eci-velocity-mag-fps </property>
<cos>
<property> flight-path/gamma-rad </property>
</cos>
</product>
<v> 2.0 </v>
</pow>
<product>
<p> guidance/executive/semimajor-axis-length-ft </p>
<p> guidance/executive/gm </p>
</product>
</quotient>
</difference>
<v> 0.5 </v>
</pow>
</function>
</fcs_function>
<fcs_function name="guidance/executive/apoapsis">
<function>
<quotient>
<difference>
<product>
<property>guidance/executive/semimajor-axis-length-ft</property>
<sum>
<value>1</value>
<property>guidance/executive/eccentricity</property>
</sum>
</product>
<property>inertial/sea-level-radius_ft</property>
</difference>
<value>5280</value>
</quotient>
</function>
</fcs_function>
<fcs_function name="guidance/executive/periapsis">
<function>
<quotient>
<difference>
<product>
<property>guidance/executive/semimajor-axis-length-ft</property>
<difference>
<value>1</value>
<property>guidance/executive/eccentricity</property>
</difference>
</product>
<property>inertial/sea-level-radius_ft</property>
</difference>
<value>5280</value>
</quotient>
</function>
</fcs_function>
<fcs_function name="guidance/executive/orbital_period">
<description>
This equation calculates orbital period in seconds. T = (2 * pi) sqrt ({a^3} / mu)
</description>)
<function>
<product>
<product>
<v> 2.0 </v>
<property> guidance/executive/pi </property>
</product>
<pow>
<quotient>
<pow>
<property>guidance/executive/semimajor-axis-length-ft</property>
<v> 3.0 </v>
</pow>
<property> guidance/executive/gm </property>
</quotient>
<v>0.5</v>
</pow>
</product>
</function>
</fcs_function>
JavaScript:
Code:
//This equation calculates the semimajor-axis length in feet. a = GM/(2*GM/r - v^2)
// GM = standard gravitational parameter of the central body, r = orbit radius (ft), v = orbit velocity (fps)
function semi_major_axis(r, v, GM){
return GM/(2.0 * GM/r - Math.pow(v,2) );
}
//This equation calculates eccentricity. e = sqrt (1 - [ r v cos(gamma) ]^2 / {a * GM})
// GM = standard gravitational parameter of the central body, r = orbit radius (ft), v = orbit velocity (fps), gamma (rad) =
function eccentricity(a, r, v, gamma, GM){
return Math.sqrt(1.0 - Math.pow( ( r * v * Math.cos(gamma) ), 2.0) / (a * GM ) );
}
//This equation calculates apoapsis in miles. ap = ((a * {1+e}) - R) / 5280
function apoapsis(a, e, R){
return ((a * (1+e))-R)/5280.0;
}
//This equation calculates periapsis in miles. pe = ((a * {1-e}) - R) / 5280
function periapsis(a, e, R){
return ((a * (1-e))-R)/5280.0;
}
//This equation calculates orbital period in seconds. T = (2 * pi) sqrt ({a^3} / GM)
function orbital_period(a, pi, GM){
return (2.0 * pi) * Math.sqrt(Math.pow(a,3.0)/GM);
}