Partial Derivatives Computation via Autodiff#
In this notebook, we show how to use the autodiff feature of \(\partial\textrm{SGP4}\). Due to the fact that it is written in pytorch
, it automatically supports automatic differentiation via torch.autograd
.
In this notebook, we show how these partial derivatives can be constructed: for more advanced examples on how to use these gradients for practical applications, see the tutorials on state_transition_matrix_computation
, covariance_propagation
, graident_based_optimization
, orbit_determination
.
import dsgp4
import torch
We create a TLE object
#as always, first, we create a TLE object:
tle=[]
tle.append('0 COSMOS 2251 DEB')
tle.append('1 34454U 93036SX 22068.91971155 .00000319 00000-0 11812-3 0 9996')
tle.append('2 34454 74.0583 280.7094 0037596 327.9100 31.9764 14.35844873683320')
tle = dsgp4.tle.TLE(tle)
print(tle)
TLE(
0 COSMOS 2251 DEB
1 34454U 93036SX 22068.91971155 .00000319 00000-0 11812-3 0 9996
2 34454 74.0583 280.7094 0037596 327.9100 31.9764 14.35844873683320
)
Now, as shown in the tle_propagation
tutorial, we can propagate the TLE. However, instead of using the standard API, we require torch.autograd
to record the operations w.r.t. the time.
Partials with respect to time#
Let’s compute the partials of the \(\partial \textrm{SGP4}\) output w.r.t. the propagation times
Single TLEs#
Let’s first see the case of single TLEs, propagated at various future times:
#let's take a random tensor of 10 tsince elements, where we track the gradients:
tsince=torch.rand((10,),requires_grad=True)
#the state is then:
state_teme = dsgp4.propagate(tle,
tsinces=tsince,
initialized=False)
#now, we can see that the gradient is tracked:
print(state_teme)
tensor([[[ 1.3832e+03, -6.9919e+03, 2.0951e+02],
[ 1.9783e+00, 5.9897e-01, 7.1995e+00]],
[[ 1.4232e+03, -6.9781e+03, 3.5617e+02],
[ 1.9468e+00, 7.5569e-01, 7.1931e+00]],
[[ 1.3264e+03, -7.0058e+03, 5.1665e+00],
[ 2.0206e+00, 3.8033e-01, 7.2029e+00]],
[[ 1.4035e+03, -6.9853e+03, 2.8367e+02],
[ 1.9625e+00, 6.7824e-01, 7.1967e+00]],
[[ 1.4314e+03, -6.9748e+03, 3.8665e+02],
[ 1.9402e+00, 7.8824e-01, 7.1914e+00]],
[[ 1.4423e+03, -6.9703e+03, 4.2699e+02],
[ 1.9313e+00, 8.3130e-01, 7.1889e+00]],
[[ 1.3561e+03, -6.9993e+03, 1.1150e+02],
[ 1.9988e+00, 4.9414e-01, 7.2019e+00]],
[[ 1.3961e+03, -6.9878e+03, 2.5675e+02],
[ 1.9683e+00, 6.4948e-01, 7.1978e+00]],
[[ 1.4217e+03, -6.9786e+03, 3.5068e+02],
[ 1.9480e+00, 7.4982e-01, 7.1934e+00]],
[[ 1.3529e+03, -7.0001e+03, 9.9830e+01],
[ 2.0013e+00, 4.8166e-01, 7.2021e+00]]],
grad_fn=<TransposeBackward0>)
Let’s now retrieve the partial derivatives of the SGP4 output w.r.t. time.
Since the state is position and velocity (i.e., \([x,y,z,v_x,v_y,v_z]\)), these partials will be all the elements of type:
Note
One thing to be careful about is that \(\partial\textrm{SGP4}\), mirroring the original \(\textrm{SGP4}\), takes the time in minutes, and returns the state in km and km/s. Hence, the derivatives will have dimensions coherent to these, and to return to SI, conversions have to be made.
partial_derivatives = torch.zeros_like(state_teme)
for i in [0,1]:
for j in [0,1,2]:
tsince.grad=None
state_teme[:,i,j].backward(torch.ones_like(tsince),retain_graph=True)
partial_derivatives[:,i,j] = tsince.grad
#let's print to screen the partials:
print(partial_derivatives)
tensor([[[ 1.1870e+02, 3.5938e+01, 4.3197e+02],
[-9.1365e-02, 4.6184e-01, -1.3875e-02]],
[[ 1.1681e+02, 4.5342e+01, 4.3159e+02],
[-9.3999e-02, 4.6089e-01, -2.3586e-02]],
[[ 1.2124e+02, 2.2820e+01, 4.3217e+02],
[-8.7625e-02, 4.6280e-01, -3.4250e-04]],
[[ 1.1775e+02, 4.0695e+01, 4.3180e+02],
[-9.2702e-02, 4.6139e-01, -1.8786e-02]],
[[ 1.1641e+02, 4.7294e+01, 4.3148e+02],
[-9.4542e-02, 4.6067e-01, -2.5604e-02]],
[[ 1.1588e+02, 4.9878e+01, 4.3133e+02],
[-9.5257e-02, 4.6036e-01, -2.8274e-02]],
[[ 1.1993e+02, 2.9649e+01, 4.3212e+02],
[-8.9581e-02, 4.6235e-01, -7.3844e-03]],
[[ 1.1810e+02, 3.8969e+01, 4.3187e+02],
[-9.2218e-02, 4.6156e-01, -1.7003e-02]],
[[ 1.1688e+02, 4.4989e+01, 4.3161e+02],
[-9.3901e-02, 4.6093e-01, -2.3222e-02]],
[[ 1.2008e+02, 2.8900e+01, 4.3213e+02],
[-8.9368e-02, 4.6241e-01, -6.6119e-03]]])
Batch TLEs#
Let’s now see how it works for batch TLEs. The API is basically identical:
#we load 6 TLEs:
inp_file="""0 PSLV DEB
1 35350U 01049QJ 22068.76869562 .00000911 00000-0 24939-3 0 9998
2 35350 98.6033 64.7516 0074531 99.8340 261.1278 14.48029442457561
0 PSLV DEB *
1 35351U 01049QK 22066.70636923 .00002156 00000-0 63479-3 0 9999
2 35351 98.8179 29.5651 0005211 45.5944 314.5671 14.44732274457505
0 SL-18 DEB
1 35354U 93014BD 22068.76520028 .00021929 00000-0 20751-2 0 9995
2 35354 75.7302 100.7819 0059525 350.7978 9.2117 14.92216400847487
0 SL-18 DEB
1 35359U 93014BJ 22068.55187275 .00025514 00000-0 24908-2 0 9992
2 35359 75.7369 156.1582 0054843 50.5279 310.0745 14.91164684775759
0 SL-18 DEB
1 35360U 93014BK 22068.44021735 .00019061 00000-0 20292-2 0 9992
2 35360 75.7343 127.2487 0071107 32.5913 327.9635 14.86997880798827
0 METEOR 2-17 DEB
1 35364U 88005Y 22067.81503681 .00001147 00000-0 84240-3 0 9995
2 35364 82.5500 92.4124 0018834 303.2489 178.0638 13.94853833332534"""
lines=inp_file.splitlines()
#let's create the TLE objects
tles=[]
for i in range(0,len(lines),3):
data=[]
data.append(lines[i])
data.append(lines[i+1])
data.append(lines[i+2])
tles.append(dsgp4.tle.TLE(data))
#we also create 9 random times, tracking the gradients:
tsinces=torch.rand((6,),requires_grad=True)
#let's propagate the batch:
state_teme = dsgp4.propagate_batch(tles,
tsinces=tsinces,
initialized=False)
Now, let’s retrieve the partial of each TLE, at each propagated time, and store them into a Nx2x3 matrix:
#let's retrieve the partials w.r.t. time:
partial_derivatives = torch.zeros_like(state_teme)
for i in [0,1]:
for j in [0,1,2]:
tsinces.grad=None
state_teme[:,i,j].backward(torch.ones_like(tsinces),retain_graph=True)
partial_derivatives[:,i,j] = tsinces.grad
#let's print to screen the partials:
print(partial_derivatives)
tensor([[[ 5.6307e+01, -3.7847e+01, 4.4360e+02],
[-2.0234e-01, -4.2667e-01, -6.7130e-03]],
[[ 1.9117e+01, -6.8198e+01, 4.4348e+02],
[-4.1179e-01, -2.3056e-01, -1.7068e-02]],
[[-1.0838e+02, -3.1776e+01, 4.4227e+02],
[ 9.6129e-02, -4.8904e-01, -1.1629e-02]],
[[-2.4392e+01, -1.1165e+02, 4.4081e+02],
[ 4.5537e-01, -1.9523e-01, -2.1684e-02]],
[[-7.7550e+01, -8.3406e+01, 4.4163e+02],
[ 3.0402e-01, -3.9214e-01, -1.8345e-02]],
[[ 4.6397e+01, -3.7284e+02, -2.3464e+02],
[ 3.8942e-02, 2.4133e-01, -3.7592e-01]]])
Partials with respect to TLE parameters#
Let’s now tackle the case in which we are interested in the partials of the \(\partial\textrm{SGP4}\) output w.r.t. the TLE parameters.
Single TLEs#
We first tackle the case of single TLE, propagated at multiple times:
In this case, we want the Jacobian of the output state, w.r.t. the following TLE parameters \(\textrm{TLE}=[n,e,i,\Omega,\omega,M,B^*,\dot{n},\ddot{n}]\), where:
\(n\) is the mean motion (also known as
no_kozai
in the original implementation) [rad/minute];\(e\) is the eccentricity [-];
\(i\) is the inclination [rad];
\(\Omega\) is the right ascension of the ascending node [rad];
\(\omega\) is the argument of perigee [rad];
\(M\) is the mean anomaly [rad];
\(B^*\) is the Bstar parameter [1/earth radii]
\(\dot{n}\) mean motion first derivative [radians/\(\textrm{minute}^2\)]
\(\ddot{n}\) mean motion second derivative [radians/\(\textrm{minute}^2\)]
#as always, first, we create a TLE object:
tle=[]
tle.append('0 COSMOS 2251 DEB')
tle.append('1 34454U 93036SX 22068.91971155 .00000319 00000-0 11812-3 0 9996')
tle.append('2 34454 74.0583 280.7094 0037596 327.9100 31.9764 14.35844873683320')
tle = dsgp4.tle.TLE(tle)
print(tle)
tle_elements=dsgp4.initialize_tle(tle,with_grad=True)
TLE(
0 COSMOS 2251 DEB
1 34454U 93036SX 22068.91971155 .00000319 00000-0 11812-3 0 9996
2 34454 74.0583 280.7094 0037596 327.9100 31.9764 14.35844873683320
)
#let's select 10 random times:
tsince=torch.rand((10,))
#and let's propagate:
state_teme=dsgp4.propagate(tle,tsince)
#now we can build the partial derivatives matrix, of shape Nx6x9 (N is the number of tsince elements, 6 is the number of elements in the state vector, and 9 is the number of elements in the TLE):
partial_derivatives = torch.zeros((len(tsince),6,9))
for k in range(len(tsince)):
for i in range(6):
tle_elements.grad=None
state_teme[k].flatten()[i].backward(retain_graph=True)
partial_derivatives[k,i,:] = tle_elements.grad
#let's print them to screen:
print(partial_derivatives)
tensor([[[-1.3116e-04, 0.0000e+00, 0.0000e+00, 9.6072e+02, 1.8732e+03,
-2.5482e+02, 1.8878e+03, -1.3765e+04, 6.9876e+03],
[-6.0700e-04, 0.0000e+00, 0.0000e+00, 6.5033e+03, 6.3590e+02,
-5.1588e+01, 6.2499e+02, 7.4664e+04, 1.3967e+03],
[-8.4522e-04, 0.0000e+00, 0.0000e+00, 7.5304e+03, 6.8527e+03,
6.9808e+01, 6.8960e+03, 1.2291e+03, 0.0000e+00],
[-1.1811e-06, 0.0000e+00, 0.0000e+00, 7.9505e-01, -1.4646e+00,
-7.0689e+00, -1.4734e+00, 9.5803e+00, -6.5182e-01],
[-7.7149e-07, 0.0000e+00, 0.0000e+00, 4.7052e+00, 7.3491e+00,
-1.3393e+00, 7.3706e+00, 7.8031e+00, 1.9678e+00],
[-4.5574e-06, 0.0000e+00, 0.0000e+00, 5.7993e+00, -2.5775e-01,
2.0504e+00, -2.7382e-01, 3.8154e+01, 0.0000e+00]],
[[-8.0466e-05, 0.0000e+00, 0.0000e+00, 9.5054e+02, 1.8910e+03,
-1.6810e+02, 1.9057e+03, -1.3885e+04, 6.9950e+03],
[-3.9600e-04, 0.0000e+00, 0.0000e+00, 6.4465e+03, 5.4572e+02,
-3.5157e+01, 5.3454e+02, 7.4580e+04, 1.3725e+03],
[-5.3393e-04, 0.0000e+00, 0.0000e+00, 7.4585e+03, 6.8553e+03,
4.4653e+01, 6.8988e+03, 7.6046e+02, 0.0000e+00],
[-7.8606e-07, 0.0000e+00, 0.0000e+00, 8.5466e-01, -1.4392e+00,
-7.0718e+00, -1.4478e+00, 9.9905e+00, -5.5741e-01],
[-4.8109e-07, 0.0000e+00, 0.0000e+00, 4.5568e+00, 7.3568e+00,
-1.3398e+00, 7.3788e+00, 5.7968e+00, 1.9865e+00],
[-3.0120e-06, 0.0000e+00, 0.0000e+00, 5.9080e+00, -1.6496e-01,
2.0512e+00, -1.8045e-01, 3.8261e+01, 0.0000e+00]],
[[-1.1341e-04, 0.0000e+00, 0.0000e+00, 9.5737e+02, 1.8792e+03,
-2.2561e+02, 1.8939e+03, -1.3804e+04, 6.9902e+03],
[-5.3543e-04, 0.0000e+00, 0.0000e+00, 6.4840e+03, 6.0553e+02,
-4.6054e+01, 5.9454e+02, 7.4633e+04, 1.3886e+03],
[-7.3773e-04, 0.0000e+00, 0.0000e+00, 7.5063e+03, 6.8537e+03,
6.1337e+01, 6.8971e+03, 1.0714e+03, 0.0000e+00],
[-1.0490e-06, 0.0000e+00, 0.0000e+00, 8.1519e-01, -1.4561e+00,
-7.0700e+00, -1.4648e+00, 9.7198e+00, -6.2003e-01],
[-6.7058e-07, 0.0000e+00, 0.0000e+00, 4.6556e+00, 7.3518e+00,
-1.3395e+00, 7.3735e+00, 7.1278e+00, 1.9741e+00],
[-4.0380e-06, 0.0000e+00, 0.0000e+00, 5.8363e+00, -2.2650e-01,
2.0507e+00, -2.4237e-01, 3.8195e+01, 0.0000e+00]],
[[-2.4109e-04, 0.0000e+00, 0.0000e+00, 9.7768e+02, 1.8393e+03,
-4.1566e+02, 1.8537e+03, -1.3555e+04, 6.9707e+03],
[-1.0112e-03, 0.0000e+00, 0.0000e+00, 6.6135e+03, 8.0302e+02,
-8.2065e+01, 7.9260e+02, 7.4884e+04, 1.4411e+03],
[-1.4854e-03, 0.0000e+00, 0.0000e+00, 7.6601e+03, 6.8449e+03,
1.1647e+02, 6.8878e+03, 2.0945e+03, 0.0000e+00],
[-1.8922e-06, 0.0000e+00, 0.0000e+00, 6.8302e-01, -1.5113e+00,
-7.0605e+00, -1.5201e+00, 8.7873e+00, -8.2673e-01],
[-1.3831e-06, 0.0000e+00, 0.0000e+00, 4.9723e+00, 7.3315e+00,
-1.3379e+00, 7.3522e+00, 1.1514e+01, 1.9323e+00],
[-7.3964e-06, 0.0000e+00, 0.0000e+00, 5.5873e+00, -4.2987e-01,
2.0481e+00, -4.4696e-01, 3.7839e+01, 0.0000e+00]],
[[-4.5902e-05, 0.0000e+00, 0.0000e+00, 9.4227e+02, 1.9043e+03,
-1.0192e+02, 1.9192e+03, -1.3979e+04, 6.9999e+03],
[-2.3787e-04, 0.0000e+00, 0.0000e+00, 6.4044e+03, 4.7685e+02,
-2.2620e+01, 4.6547e+02, 7.4533e+04, 1.3538e+03],
[-3.1239e-04, 0.0000e+00, 0.0000e+00, 7.4028e+03, 6.8565e+03,
2.5458e+01, 6.9001e+03, 4.0219e+02, 0.0000e+00],
[-4.7914e-07, 0.0000e+00, 0.0000e+00, 8.9977e-01, -1.4196e+00,
-7.0732e+00, -1.4282e+00, 1.0295e+01, -4.8533e-01],
[-2.7838e-07, 0.0000e+00, 0.0000e+00, 4.4415e+00, 7.3618e+00,
-1.3400e+00, 7.3842e+00, 4.2634e+00, 2.0006e+00],
[-1.8263e-06, 0.0000e+00, 0.0000e+00, 5.9883e+00, -9.4147e-02,
2.0516e+00, -1.0920e-01, 3.8312e+01, 0.0000e+00]],
[[-9.5022e-05, 0.0000e+00, 0.0000e+00, 9.5367e+02, 1.8857e+03,
-1.9409e+02, 1.9004e+03, -1.3848e+04, 6.9929e+03],
[-4.5877e-04, 0.0000e+00, 0.0000e+00, 6.4633e+03, 5.7275e+02,
-4.0082e+01, 5.6166e+02, 7.4603e+04, 1.3798e+03],
[-6.2474e-04, 0.0000e+00, 0.0000e+00, 7.4802e+03, 6.8546e+03,
5.2193e+01, 6.8981e+03, 9.0104e+02, 0.0000e+00],
[-9.0533e-07, 0.0000e+00, 0.0000e+00, 8.3686e-01, -1.4468e+00,
-7.0710e+00, -1.4555e+00, 9.8689e+00, -5.8572e-01],
[-5.6520e-07, 0.0000e+00, 0.0000e+00, 4.6016e+00, 7.3546e+00,
-1.3397e+00, 7.3765e+00, 6.3985e+00, 1.9809e+00],
[-3.4762e-06, 0.0000e+00, 0.0000e+00, 5.8758e+00, -1.9277e-01,
2.0510e+00, -2.0844e-01, 3.8234e+01, 0.0000e+00]],
[[-1.7153e-06, 0.0000e+00, 0.0000e+00, 9.2938e+02, 1.9236e+03,
-4.7098e+00, 1.9386e+03, -1.4124e+04, 7.0058e+03],
[-9.6720e-06, 0.0000e+00, 0.0000e+00, 6.3445e+03, 3.7564e+02,
-4.2038e+00, 3.6395e+02, 7.4490e+04, 1.3262e+03],
[-1.2185e-05, 0.0000e+00, 0.0000e+00, 7.3197e+03, 6.8571e+03,
-2.7388e+00, 6.9009e+03, -1.2457e+02, 0.0000e+00],
[-1.9888e-08, 0.0000e+00, 0.0000e+00, 9.6538e-01, -1.3906e+00,
-7.0740e+00, -1.3991e+00, 1.0730e+01, -3.7938e-01],
[-1.0659e-08, 0.0000e+00, 0.0000e+00, 4.2691e+00, 7.3680e+00,
-1.3401e+00, 7.3908e+00, 2.0079e+00, 2.0208e+00],
[-7.5225e-08, 0.0000e+00, 0.0000e+00, 6.1021e+00, 9.8687e-03,
2.0517e+00, -4.5251e-03, 3.8340e+01, 0.0000e+00]],
[[-2.0237e-05, 0.0000e+00, 0.0000e+00, 9.3518e+02, 1.9152e+03,
-4.7599e+01, 1.9301e+03, -1.4059e+04, 7.0034e+03],
[-1.0978e-04, 0.0000e+00, 0.0000e+00, 6.3706e+03, 4.2031e+02,
-1.2329e+01, 4.0875e+02, 7.4505e+04, 1.3384e+03],
[-1.4093e-04, 0.0000e+00, 0.0000e+00, 7.3565e+03, 6.8570e+03,
9.7014e+00, 6.9007e+03, 1.0787e+02, 0.0000e+00],
[-2.2373e-07, 0.0000e+00, 0.0000e+00, 9.3653e-01, -1.4034e+00,
-7.0738e+00, -1.4120e+00, 1.0540e+01, -4.2613e-01],
[-1.2434e-07, 0.0000e+00, 0.0000e+00, 4.3456e+00, 7.3655e+00,
-1.3401e+00, 7.3881e+00, 3.0034e+00, 2.0119e+00],
[-8.4912e-07, 0.0000e+00, 0.0000e+00, 6.0525e+00, -3.6023e-02,
2.0517e+00, -5.0709e-02, 3.8335e+01, 0.0000e+00]],
[[-8.0360e-05, 0.0000e+00, 0.0000e+00, 9.5051e+02, 1.8910e+03,
-1.6790e+02, 1.9058e+03, -1.3885e+04, 6.9950e+03],
[-3.9553e-04, 0.0000e+00, 0.0000e+00, 6.4464e+03, 5.4551e+02,
-3.5121e+01, 5.3434e+02, 7.4580e+04, 1.3724e+03],
[-5.3326e-04, 0.0000e+00, 0.0000e+00, 7.4584e+03, 6.8553e+03,
4.4597e+01, 6.8988e+03, 7.5941e+02, 0.0000e+00],
[-7.8517e-07, 0.0000e+00, 0.0000e+00, 8.5480e-01, -1.4391e+00,
-7.0718e+00, -1.4478e+00, 9.9914e+00, -5.5720e-01],
[-4.8048e-07, 0.0000e+00, 0.0000e+00, 4.5565e+00, 7.3568e+00,
-1.3398e+00, 7.3788e+00, 5.7923e+00, 1.9866e+00],
[-3.0085e-06, 0.0000e+00, 0.0000e+00, 5.9082e+00, -1.6475e-01,
2.0512e+00, -1.8024e-01, 3.8261e+01, 0.0000e+00]],
[[-2.1478e-04, 0.0000e+00, 0.0000e+00, 9.7411e+02, 1.8469e+03,
-3.7974e+02, 1.8614e+03, -1.3601e+04, 6.9749e+03],
[-9.1935e-04, 0.0000e+00, 0.0000e+00, 6.5884e+03, 7.6572e+02,
-7.5259e+01, 7.5519e+02, 7.4827e+04, 1.4313e+03],
[-1.3353e-03, 0.0000e+00, 0.0000e+00, 7.6316e+03, 6.8469e+03,
1.0605e+02, 6.8900e+03, 1.9018e+03, 0.0000e+00],
[-1.7359e-06, 0.0000e+00, 0.0000e+00, 7.0820e-01, -1.5009e+00,
-7.0627e+00, -1.5098e+00, 8.9680e+00, -7.8769e-01],
[-1.2384e-06, 0.0000e+00, 0.0000e+00, 4.9136e+00, 7.3358e+00,
-1.3383e+00, 7.3567e+00, 1.0686e+01, 1.9403e+00],
[-6.7656e-06, 0.0000e+00, 0.0000e+00, 5.6358e+00, -3.9143e-01,
2.0487e+00, -4.0830e-01, 3.7923e+01, 0.0000e+00]]])
Batch TLEs:#
As for the time derivatives, the API stays practically identical:
#we load 6 TLEs:
inp_file="""0 PSLV DEB
1 35350U 01049QJ 22068.76869562 .00000911 00000-0 24939-3 0 9998
2 35350 98.6033 64.7516 0074531 99.8340 261.1278 14.48029442457561
0 PSLV DEB *
1 35351U 01049QK 22066.70636923 .00002156 00000-0 63479-3 0 9999
2 35351 98.8179 29.5651 0005211 45.5944 314.5671 14.44732274457505
0 SL-18 DEB
1 35354U 93014BD 22068.76520028 .00021929 00000-0 20751-2 0 9995
2 35354 75.7302 100.7819 0059525 350.7978 9.2117 14.92216400847487
0 SL-18 DEB
1 35359U 93014BJ 22068.55187275 .00025514 00000-0 24908-2 0 9992
2 35359 75.7369 156.1582 0054843 50.5279 310.0745 14.91164684775759
0 SL-18 DEB
1 35360U 93014BK 22068.44021735 .00019061 00000-0 20292-2 0 9992
2 35360 75.7343 127.2487 0071107 32.5913 327.9635 14.86997880798827
0 METEOR 2-17 DEB
1 35364U 88005Y 22067.81503681 .00001147 00000-0 84240-3 0 9995
2 35364 82.5500 92.4124 0018834 303.2489 178.0638 13.94853833332534"""
lines=inp_file.splitlines()
#let's create the TLE objects
tles=[]
for i in range(0,len(lines),3):
data=[]
data.append(lines[i])
data.append(lines[i+1])
data.append(lines[i+2])
tles.append(dsgp4.tle.TLE(data))
#we also create 9 random times, tracking the gradients:
tsinces=torch.rand((6,))
#let's now initialize the TLEs, activating the gradient tracking for the TLE parameters:
tle_elements=dsgp4.initialize_tle(tles,with_grad=True)
#let's now propagate the batch of TLEs:
state_teme = dsgp4.propagate_batch(tles,tsinces)
Finally, we can build the matrix that contains the partial of the SGP4 output w.r.t. the TLE parameters, for each TLE:
#now we can build the partial derivatives matrix, of shape Nx6x9 (N is the number of tsince elements, 6 is the number of elements in the state vector, and 9 is the number of elements in the TLE):
partial_derivatives = torch.zeros((len(tsinces),6,9))
for k in range(len(tsinces)):
for i in range(6):
tle_elements[k].grad=None
state_teme[k].flatten()[i].backward(retain_graph=True)
partial_derivatives[k,i,:] = tle_elements[k].grad
#let's print them to screen:
print(partial_derivatives)
tensor([[[-2.8317e-04, 0.0000e+00, 0.0000e+00, -1.2881e+03, 8.6195e+02,
2.0582e+02, 8.3764e+02, -3.1920e+04, -6.4215e+03],
[-1.3691e-03, 0.0000e+00, 0.0000e+00, 2.2039e+03, -6.7015e+02,
-9.9142e+01, -7.1619e+02, -6.8056e+04, 3.0665e+03],
[ 2.1658e-03, 0.0000e+00, 0.0000e+00, -1.3905e+04, 7.0376e+03,
-3.2203e+01, 7.0221e+03, 1.0482e+03, 0.0000e+00],
[-9.6130e-07, 0.0000e+00, 0.0000e+00, 3.0804e+00, -3.2321e+00,
6.6798e+00, -3.2224e+00, 3.0085e+00, 7.5314e-01],
[-8.3334e-07, 0.0000e+00, 0.0000e+00, 6.7890e+00, -6.7489e+00,
-3.1516e+00, -6.7477e+00, -7.3652e+00, 8.8022e-01],
[-3.3875e-06, 0.0000e+00, 0.0000e+00, -7.2485e-01, -2.9533e-01,
-1.1149e+00, -2.4047e-01, 3.8877e+01, 0.0000e+00]],
[[-2.5369e-04, 0.0000e+00, 0.0000e+00, -4.9544e+03, 3.3350e+02,
1.0982e+02, 3.3152e+02, -6.5440e+04, -3.4814e+03],
[-2.7291e-04, 0.0000e+00, 0.0000e+00, -1.0472e+03, -1.0663e+03,
-1.9748e+02, -1.0684e+03, -3.7307e+04, 6.2079e+03],
[ 7.2307e-04, 0.0000e+00, 0.0000e+00, -9.8808e+03, 7.0349e+03,
-3.2651e+01, 7.0401e+03, 1.0378e+03, 0.0000e+00],
[-6.0837e-07, 0.0000e+00, 0.0000e+00, 4.7661e+00, -6.5324e+00,
3.6442e+00, -6.5347e+00, -1.4230e+00, 1.1205e+00],
[ 3.4086e-07, 0.0000e+00, 0.0000e+00, 1.7212e+00, -3.6628e+00,
-6.4268e+00, -3.6646e+00, -7.7365e+00, 3.4737e-01],
[-3.8379e-06, 0.0000e+00, 0.0000e+00, 5.5126e+00, -2.4079e-01,
-1.1434e+00, -2.3822e-01, 3.8979e+01, 0.0000e+00]],
[[-1.2864e-03, 0.0000e+00, 0.0000e+00, 6.6485e+02, -1.6011e+03,
3.7907e+02, -1.6217e+03, 1.2854e+04, -6.7791e+03],
[ 7.3426e-03, 0.0000e+00, 0.0000e+00, -6.9669e+03, -7.1789e+02,
7.5376e+01, -7.1743e+02, -6.9970e+04, -1.3908e+03],
[-4.2920e-04, 0.0000e+00, 0.0000e+00, 2.5542e+03, 6.7076e+03,
9.4244e+01, 6.7866e+03, 1.8256e+03, 0.0000e+00],
[ 6.9113e-06, 0.0000e+00, 0.0000e+00, -1.3789e+00, 1.5261e+00,
7.2272e+00, 1.5373e+00, -7.6464e+00, 7.7765e-01],
[ 2.6250e-06, 0.0000e+00, 0.0000e+00, -2.3984e+00, -7.4501e+00,
1.3787e+00, -7.4928e+00, -1.0436e+01, -1.7565e+00],
[-2.8588e-05, 0.0000e+00, 0.0000e+00, 7.0947e+00, -4.1552e-01,
1.8669e+00, -4.2749e-01, 3.7336e+01, 0.0000e+00]],
[[-1.1855e-03, 0.0000e+00, 0.0000e+00, 5.0368e+03, -5.7998e+02,
4.2936e+01, -5.5752e+02, 6.5065e+04, -2.7828e+03],
[-2.2255e-04, 0.0000e+00, 0.0000e+00, 6.4789e+02, -1.6171e+03,
1.0470e+02, -1.6404e+03, -2.8873e+04, -6.3677e+03],
[ 2.6842e-03, 0.0000e+00, 0.0000e+00, -1.0333e+04, 6.7362e+03,
2.4874e+01, 6.7841e+03, 4.3109e+02, 0.0000e+00],
[ 2.0304e-06, 0.0000e+00, 0.0000e+00, -5.7059e+00, 6.9628e+00,
2.9689e+00, 6.9853e+00, -1.3786e+00, 1.7757e+00],
[ 1.0542e-06, 0.0000e+00, 0.0000e+00, 1.1628e+00, -3.0344e+00,
6.7233e+00, -3.0527e+00, -9.8214e+00, -6.0325e-01],
[-7.0105e-06, 0.0000e+00, 0.0000e+00, 4.8925e+00, -1.5360e-01,
1.8638e+00, -1.2348e-01, 3.7658e+01, 0.0000e+00]],
[[-1.1004e-03, 0.0000e+00, 0.0000e+00, 4.9923e+03, -1.3088e+03,
6.2978e+01, -1.3089e+03, 4.3118e+04, -5.5162e+03],
[ 8.1120e-04, 0.0000e+00, 0.0000e+00, -3.5440e+03, -1.1088e+03,
5.1693e+01, -1.1431e+03, -5.6818e+04, -4.2201e+03],
[ 1.5130e-03, 0.0000e+00, 0.0000e+00, -7.1883e+03, 6.7324e+03,
1.6832e+01, 6.8142e+03, 2.6956e+02, 0.0000e+00],
[ 1.5522e-06, 0.0000e+00, 0.0000e+00, -3.6545e+00, 4.6281e+00,
5.8587e+00, 4.6510e+00, -6.4417e+00, 1.2340e+00],
[ 1.8008e-07, 0.0000e+00, 0.0000e+00, 2.1535e+00, -6.0379e+00,
4.4575e+00, -6.0794e+00, -7.3697e+00, -1.4127e+00],
[-5.2812e-06, 0.0000e+00, 0.0000e+00, 6.3193e+00, -1.1594e-01,
1.8672e+00, -8.9007e-02, 3.7839e+01, 0.0000e+00]],
[[-2.1852e-06, 0.0000e+00, 0.0000e+00, -6.0414e+02, 7.6226e+02,
6.1364e+03, 7.5937e+02, 7.1406e+03, 3.8957e+03],
[ 7.0910e-06, 0.0000e+00, 0.0000e+00, -4.1733e+03, -6.1774e+03,
2.5774e+02, -6.1543e+03, 4.1569e+04, -6.3974e+02],
[ 1.3204e-05, 0.0000e+00, 0.0000e+00, 5.9555e+03, -3.8334e+03,
8.0336e+02, -3.8187e+03, -6.7981e+04, 0.0000e+00],
[-6.1152e-08, 0.0000e+00, 0.0000e+00, -7.5487e-01, 6.4732e-01,
-3.8657e+00, 6.4607e-01, 4.3129e+00, 6.2391e+00],
[ 4.9132e-07, 0.0000e+00, 0.0000e+00, 6.3317e+00, 3.9354e+00,
-1.6148e-01, 3.9283e+00, -3.3484e+01, 7.6920e-01],
[ 3.0852e-07, 0.0000e+00, 0.0000e+00, 3.7317e+00, -6.2164e+00,
-5.0426e-01, -6.2046e+00, -2.2245e+01, 0.0000e+00]]])