Field propagation not returning a DerivativeStructure


#1

Hi all!

I am attempting to obtain the partial derivatives of an orbit propagation without much success. The propagation works fine, when retrieving the final position, instead of a DerivativeStructure I get a RealFieldElement which does not have (or I don’t know how to access) the getAllDerivatives() method.

I am using the python wrapper.

Here I attach a minimalist example where I generate a FieldKeplerianOrbit but when retrieving the mean motion of this orbit I get a RealFieldElement instead of DerivativeStructure. A simple hand calculation provides the desired DerivativeStructure class. Is there anything that I am missing?

import orekit
from org.orekit.utils import Constants
from org.orekit.frames import FramesFactory
from org.orekit.time import FieldAbsoluteDate
from org.orekit.orbits import FieldKeplerianOrbit
from org.orekit.orbits import PositionAngle
from org.hipparchus.analysis.differentiation import DSFactory
from math import radians

# Initialize JVM
vm = orekit.initVM()

# DS Factory
factory = DSFactory(2, 2)

# Some Constants
ae = Constants.WGS84_EARTH_EQUATORIAL_RADIUS
mu = Constants.WGS84_EARTH_MU

# Orbit parameters
ra = 500 * 1000.0         # km Apogee!
rp = 815 * 1000.0         # km Perigee!
i_nominal = radians(98.55)      # inclination
omega_nominal = radians(90.0)   # perigee argument
raan_nominal = radians(5.1917)  # right ascension of ascending node
lv_nominal = radians(0.0567634)    # True anomaly
a_nominal = (rp + ra + 2 * ae) / 2.0    # semi major axis in km
e_nominal = 1.0 - (rp + ae) / a_nominal  # Eccentricity

# Generate field
a = factory.variable(0, a_nominal)
i = factory.variable(1, i_nominal)
omega = factory.constant(omega_nominal)
raan = factory.constant(raan_nominal)
lv = factory.constant(lv_nominal)
e = factory.constant(e_nominal)
field = a.getField()
zero = field.getZero()

epochDate = FieldAbsoluteDate(field)

# Inertial frame where the satellite is defined
inertialFrame = FramesFactory.getEME2000()

# Orbit construction as Keplerian
orbit = FieldKeplerianOrbit(a, e, i, omega, raan, lv,
                            PositionAngle.TRUE,
                            inertialFrame, epochDate, mu)

# Calulate mean motion
n_handcalc = a.abs().reciprocal().multiply(mu).sqrt().divide(a.abs())
# Print class
print n_handcalc.class_
# output: class org.hipparchus.analysis.differentiation.DerivativeStructure

# Mean motion via FieldKeplerianOrbit
n_orbit = orbit.getKeplerianMeanMotion()
# Print class
print n_orbit.class_
# output: interface org.hipparchus.RealFieldElement


# n_handcalc has getAllDerivatives() method but n_orbit not, why?

Thanks!


#2

Solved!

Just need to add .of_(DerivativeStructure) when creating new instances.

import orekit
from org.orekit.utils import Constants
from org.orekit.frames import FramesFactory
from org.orekit.time import FieldAbsoluteDate
from org.orekit.orbits import FieldKeplerianOrbit
from org.orekit.orbits import PositionAngle
from org.hipparchus.analysis.differentiation import DSFactory
from org.hipparchus.analysis.differentiation import DerivativeStructure
from math import radians

# Initialize JVM
vm = orekit.initVM()

# DS Factory
factory = DSFactory(2, 2)

# Some Constants
ae = Constants.WGS84_EARTH_EQUATORIAL_RADIUS
mu = Constants.WGS84_EARTH_MU

# Orbit parameters
ra = 500 * 1000.0         # km Apogee!
rp = 815 * 1000.0         # km Perigee!
i_nominal = radians(98.55)      # inclination
omega_nominal = radians(90.0)   # perigee argument
raan_nominal = radians(5.1917)  # right ascension of ascending node
lv_nominal = radians(0.0567634)    # True anomaly
a_nominal = (rp + ra + 2 * ae) / 2.0    # semi major axis in km
e_nominal = 1.0 - (rp + ae) / a_nominal  # Eccentricity

# Field stuff
a = factory.variable(0, a_nominal)
i = factory.variable(1, i_nominal)
omega = factory.constant(omega_nominal)
raan = factory.constant(raan_nominal)
lv = factory.constant(lv_nominal)
e = factory.constant(e_nominal)
field = a.getField().of_(DerivativeStructure)
zero = field.getZero()

epochDate = FieldAbsoluteDate(field).of_(DerivativeStructure)

# Inertial frame where the satellite is defined
inertialFrame = FramesFactory.getEME2000()

# Orbit construction as Keplerian
orbit = FieldKeplerianOrbit(a, e, i, omega, raan, lv,
                            PositionAngle.TRUE,
                            inertialFrame, epochDate, mu).of_(DerivativeStructure)

# Calulate mean motion
n_handcalc = a.abs().reciprocal().multiply(mu).sqrt().divide(a.abs())
# Print class
print n_handcalc.class_
# output: class org.hipparchus.analysis.differentiation.DerivativeStructure

# Mean motion via FieldKeplerianOrbit
n_orbit = orbit.getKeplerianMeanMotion()
# Print classs
print n_orbit.class_
# output: class org.hipparchus.analysis.differentiation.DerivativeStructure

#3

Hi,

Great you found a way to solve it!

Another way is to cast the n_orbit to a DerivativeStructure if you know that that’s the Interface you want to have.

The syntax for this is:

In[6]: DerivativeStructure.cast_(n_orbit)
Out[5]: <DerivativeStructure: org.hipparchus.analysis.differentiation.DerivativeStructure@9878fede>

Regards
/Petrus