TimeEvaluables¶
In this notebook we are going to learn more about TimeEvaluables in GODOT.
A TimeEvaluable class is simply a class which contains an evaluation method which takes a time and returns a value.
TimeEvaluable is actually a family of interfaces:
- ScalarTimeEvaluable
- Vector3TimeEvaluable
- Vector6TimeEvaluable
- Vector7TimeEvaluable
- Vector9TimeEvaluable
There are other similar interfaces you can learn about later. This is where to get started.
There are many implementations of actual classes (not just interfaces) which implement the TimeEvaluable interfaces: look for names containing https://orbit_software.io.esa.int/godotdocs/0.3.1/api/godot.model.html#module-godot.model.common in the godot common library, or peruse the list of (C++) derived types.
This is a big part of the GODOT computation philosophy.
Here we are going to explain how TimeEvaluables can be used, and extended. This all applies to C++ as well as Python of course.
Importing Modules¶
First we import some necessary modules:
from godot.core import tempo
import godot.core.autodif as ad
import godot.model.common as common
# optionally avoid verbose logging messages
import godot.core.util as util
util.suppressLogger()
Next we will define our first ScalarTimeEvaluable class:
class Square ( common.ScalarTimeEvaluable ):
def __init__(self, x : common.ScalarTimeEvaluable ):
super().__init__()
self.__x = x
def eval(self, e ):
x = self.__x.eval(e)
return x * x
This simple takes computes the square of its input. The input must also be a ScalarTimeEvaluable.
The eval method evaluates its inputs and does the computation.
We introduce common.ConstantScalarTimeEvaluable, which is simply a constant wrapped in a constant TimeEvaluable:
x = common.ConstantScalarTimeEvaluable(2.0)
z = x.eval(tempo.Epoch())
print(z)
Now we create an instance of our Square object, and evaluate it:
f = Square(x)
f.eval(tempo.Epoch())
Now these objects are very useful.
autodif¶
For one thing, I can use the same class to compute with and without gradients. Suppose I pass Square an autodif type, then I can evaluate it with or without partials by using Epoch or XEpoch respectively.
x = common.ConstantScalarTimeEvaluable(ad.Scalar( 2.0, "X" ))
f = Square(x)
z = f.eval(tempo.Epoch())
print(z)
z = f.eval(tempo.XEpoch())
print(z)
Note that to use autodif properly you need to use the bridge module for Python.
TimeEvaluable Operators¶
For another thing, I can use some mathematical operations on these objects:
x = common.ConstantScalarTimeEvaluable(2.0)
f = Square(x)
g = f + 3.0 * f
z = g.eval(tempo.Epoch())
print(z)
z = g.eval(tempo.XEpoch())
print(z)
x = common.ConstantScalarTimeEvaluable(ad.Scalar( 2.0, "X" ))
f = Square(x)
g = f + x * f
z = g.eval(tempo.Epoch())
print(z)
z = g.eval(tempo.XEpoch())
print(z)
Now it is important to realise that all this makes use of 'lazy evaluation'. The operations do not actually do the computation, they only set up a tree of functions which can all later, evaluate one another.
To illustrate this, here is another example.
class OtherSquare ( common.ScalarTimeEvaluable ):
def __init__(self, x : common.ScalarTimeEvaluable ):
super().__init__()
self.setx(x)
def setx(self,x : common.ScalarTimeEvaluable):
self.__x = x
def eval(self, e ):
x = self.__x.eval(e)
return x * x
Lets create an OtherSquare object, initialise with value 2.0, and compute a function of f.
x = common.ConstantScalarTimeEvaluable(ad.Scalar( 2.0, "X" ))
f = OtherSquare(x)
g = f + x * f
z = g.eval(tempo.Epoch())
print(z)
z = g.eval(tempo.XEpoch())
print(z)
Now, we set the value of x inside f to be 3.0, and evaluate g. Because of the lazy evaluation, The result of g is now updated. In addition, there is a another gradient element.
f.setx(common.ConstantScalarTimeEvaluable(ad.Scalar( 3.0, "newX" )))
z = g.eval(tempo.XEpoch())
print(z)
The reason there are now 2 partials is because g = f + x * f, there is the gradient of the X in that expression, and the gradient of the newX which is introduced to f.
So, it is very powerful but you have to be careful. When generating you own TimeEvaluable classes for use with autodif, always write unit tests to compute the gradients numerically.
Gradients with respect to initial state¶
The following is an interesting example of how to test the partials, using the functionality just described above.
First, we load cosmos:
from godot import cosmos
from godot.core import util
Next, we load the universe configuration and create a universe object:
# create the universe
uniConfig = cosmos.util.load_yaml('universe.yml')
uni = cosmos.Universe(uniConfig)
Now we create the trajectory:
# create the trajectory
traConfig = cosmos.util.load_yaml('trajectory.yml')
tra = cosmos.Trajectory(uni, traConfig)
Now, propagate the trajectory with all available partial derivatives tracked:
# avoid verbose propagator logging
util.suppressLogger()
# track all available partials
uni.parameters.track(uni.parameters.contains())
print(uni.parameters)
# propagate with partials
tra.compute(partials = True)
We now create a TimeEValuable class to compute the distance of the propagated spacecraft to Earth
# define a function
class f(common.ScalarTimeEvaluable):
def __init__(self):
super().__init__()
pass
def eval(self, e):
return uni.frames.distance('Earth', 'GeoSat_center', e)
Now, we evaluate it and compute the gradient with respect to the initial state time.
# print out the computed partials with respect to the initial epoch, "launch_dt"
func = f()
e = tempo.XEpoch("2000-01-01T12:00:00.000 TDB")
d = common.resolve( func.eval(e) )
df_dt_ad = d.at("launch_dt")[0]
print(df_dt_ad)
I want to test the value here using numerical differentiation.
# set the initial parameter value to -0.5 secs, and evaluate the function at epoch
uni.parameters.get("launch_dt").setPhysicalValue(-0.5)
tra.compute(partials = False)
fminus = func.eval( e.value() )
# set the initial parameter value to +0.5 secs, and evaluate the function at epoch
uni.parameters.get("launch_dt").setPhysicalValue(+0.5)
tra.compute(partials = False)
fplus = func.eval( e.value() )
# compute the gradient numerically
df_dt_num = fplus-fminus
# print the difference
print( 'Autodif: %15.6f\nNumerical: %15.6f\nDiff(percent): %15.6f' % ( df_dt_ad, df_dt_num, 100 * ( df_dt_ad - df_dt_num )/df_dt_num ) )