Information about DSST propagation tutorial

Orekit usage
1

Dear all,

I was running the DSSTPropagation tutorial, and I got some unexpected results (at least from my perspective) when trying to compare the propagated state of the DSST propagator with the corresponding one from the numerical propagator.

I am using the code in the DSSTPropagation class as it is, without applying any modification to it.
The only changes I had to implement were in the associated .yaml file holding the inputs for the simulation. Specifically, I’m setting in the .yaml both “initialOrbitIsOsculating” and “outputOrbitIsOsculating” to true, to enable a meaningful comparison with respect to numerical propagation output. Morevoer, I run few test cases, using in each a different setup of the force models used (see later for details).

I am providing below the results for each test case, where for convenience I am reporting only the position and velocity in km at a given date I selected for testing (14 days from propagation start epoch).

Test case 1 (only 6x6 geopotential, corresponding to default setup):

DSST

time from start (s) position along X (km) position along Y (km) position along Z (km) velocity along X (km) velocity along Y (km) velocity along Z (km)
1209600.00000000000 -2.5453788793846297e+03 -5.1763794489020730e+03 4.3013005597937070e+03 -3.1987758999384988e+00 -3.3047296177848800e+00 -5.8533803459324270e+00

Numerical

time from start (s) position along X (km) position along Y (km) position along Z (km) velocity along X (km) velocity along Y (km) velocity along Z (km)
1209600.00000000000 -2.5396764434146410e+03 -5.1600673432357220e+03 4.3241780097195015e+03 -3.2151893494493330e+00 -3.3212013494010140e+00 -5.8350406192403620e+00

DSST execution time (without large file write) : 2.719

Numerical execution time (including output): 1.553

Difference in position (norm of delta position vector)
28.67 km

Difference in velocity (norm of delta velocity vector)
0.029615 km/s

Test case 2 (full force model with 40x40 geopotential, Sun+Moon third body, drag (area=5.0 m^2, cd=2.0), SRP (area=5.0 m^2, cr=1.2)):

DSST

time from start (s) position along X (km) position along Y (km) position along Z (km) velocity along X (km) velocity along Y (km) velocity along Z (km)
1209600.0000000000 -2.5500687511188790e+03 -5.1795792955146040e+03 4.2943581108833205e+03 -3.1960697147360935e+00 -3.2978478046877710e+00 -5.8589912320143800e+00

Numerical

time from start (s) position along X (km) position along Y (km) position along Z (km) velocity along X (km) velocity along Y (km) velocity along Z (km)
1209600.00000000000 -2.4529313875578036e+03 -5.0678381035775380e+03 4.4795981857412800e+03 -3.2875023751665213e+00 -3.4678823032998585e+00 -5.7086061329577660e+00

DSST execution time (without large file write) : 31.607

Numerical execution time (including output): 23.363

Difference in position (norm of delta position vector)
237.14 km

Difference in velocity (norm of delta velocity vector)
0.2447 km/s

Test case 3 (force model with Sun+Moon third body, drag (area=5.0 m^2, cd=2.0), SRP (area=5.0 m^2, cr=1.2)):

DSST

time from start (s) position along X (km) position along Y (km) position along Z (km) velocity along X (km) velocity along Y (km) velocity along Z (km)
1209600.000000000 5.3884825614808080e+02 -9.9390158981442240e+02 7.1067305889496430e+03 -5.4441994291983740e+00 -5.0721742075551990e+00 - 2.9957683042493455e-01

Numerical

time from start (s) position along X (km) position along Y (km) position along Z (km) velocity along X (km) velocity along Y (km) velocity along Z (km)
1209600.000000000 5.3874220520087060e+02 -9.9400037501156190e+02 7.1067246353801500e+03 -5.4442107504746510e+00 -5.0721534238986300e+00 -2.9972492472615525e-01

DSST execution time (without large file write) : 3.518

Numerical execution time (including output): 15.961

Difference in position (norm of delta position vector)
0.1451 km

Difference in velocity (norm of delta velocity vector)
1.4997e-04 km/s

As you can see above, I am considering also a test case where using all the force models that can be set in DSST propagator, but neglecting the geopotential, as I could see that there are few bug issues still open about DSST propagator, where the geopotential force model seems to be affected (so to avoid considering in results the possible effect of bugs still open).
Clearly, all the force models are set similarly in DSST and numerical propagator before making any comparison.

As you can see from the results above, after only 14 days of propagation, in the first two test cases I get errors in propagated position of several kilometers (increasing to hundreds on kilometers in second test case), and errors in propagated velocity of the order of fraction of km/s, but still increasing in second test case.
Are these errors expected (at least in terms of order of magnitude with respect to numerical propagation results)?
I would have expected much smaller differences, but maybe I am wrong.

However, in third test case these errors are way much smaller. Therefore I am suspecting that the discrepancies may be due to the bugs in the DSST geopotential model, can it be?

Moreover, in first two test cases (as also reported in results above) I could see DSST computational times higher with respect to numerical propagator computational times, wheres I would have expected the opposite, i.e. DSST to be faster than numerical to some extent.
Only in third test case I can see the expected behaviour in terms of computational times, i.e. DSST being almost x5 times faster than numerical.

Can this behaviour be also due to the open bugs in DSST proapgator?

Again, just to remark that I did not make any change in the tutorial code (i.e. in the Java code), but I changed only few settings in the associated .yaml configuration file, as reported before. Therefore, I am assuming that all the tutorial logic is implemented correctly, and both DSST and numerical propagators are set up correctly.

Looking forward for your feedback about these behaviours, I would like to thank you in advance for your time!

2

Hi there,

@MaximeJ improved Orekit’s DSST accuracy (on zonal harmonics) in version 12.1 or thereabout, but the tutorials has not been modified to call the improved stuff.

Try this:
dsstProp.addForceModel(new DSSTZonal(earthFrame, unnormalized));

About performance, well I think there is room for improvement in Orekit’s DSST.
The harmonics seem to be a bit slow indeed, probably due to a heavy use of automatic differentiation (as opposed to analytical formulas).

Cheers,
Romain.

3

Hi,

Thank you @Serrof for your first feedback!

@nick, I suggest you use a TOD here as Earth frame since it has precession and nutation but is faster to compute than ITRF. It’s generally a good performance/precision balance.
It should reduce the discrepancy in inclination and RAAN, but not the one in SMA and AOL I’m afraid…

I think so, yes, unfortunately… . I suspect #1491 is the main culprit, with the J2 zonal harmonic being the biggest contributor.
The test case that led to the discovery of this issue is the same type of orbit (LEO SSO frozen orbit).
Could you try with only J2 and check the differences?
Also, could you please plot the difference in circular elements so that we can make sure that the error is consistent with what is reported in #1491? I may have a script somewhere to do it if you need.

Here I think this is due to the tesseral harmonics. Maybe you can play around with the configuration but it’s a bit hard to know exactly what to do. And as long as #1491 is there you won’t be able to see the difference in precision since J2 contribution is much more important than the tesserals.

Cheers,
Maxime

4

Hi all,

thanks a lot for sharing your feedback, really appreciated!

I tried to setup the DSST zonal force model using

dsstProp.addForceModel(new DSSTZonal(earthFrame, unnormalized));

with earthFrame = FramesFactory.getTOD(IERSConventions.IERS_2010, false).

@MaximeJ, setting up the propagation with the change above and simulating only a 2x0 gravity field, I get the following results:

DSST

time from start (s) position along X (km) position along Y (km) position along Z (km) velocity along X (km) velocity along Y (km) velocity along Z (km)
1209600.0000000000 -2.5726574696566240e+03 -5.2001561516874890e+03 4.2650601665415030e+03 -3.1826275665848470e+00 -3.2630152975253353e+00 -5.8786508192580040e+00

Numerical

time from start (s) position along X (km) position along Y (km) position along Z (km) velocity along X (km) velocity along Y (km) velocity along Z (km)
1209600.0000000000 -2.5622512743220530e+03 -5.1880208524541520e+03 4.2860176210133960e+03 -3.1929963020466480e+00 -3.2821900672145605e+00 -5.8623346517843800e+00

DSST execution time (without large file write) : 1.643

Numerical execution time (including output): 1.564

Difference in position (norm of delta position vector)
26.358 km

Difference in velocity (norm of delta velocity vector)
0.027229 km/s

As you can notice, the discrepancy from the 6x6 geopotential test case (in my previous post) is decreased a bit, but still quite similar.

(@MaximeJ, I am sorry, but unfortunately I am not able to provide you the plots right now, as limited by my current machine setup)

I hope that following future fixes of the open bugs on DSST propagator, we can have this propagator as full competitor of the numerical in applications where computational speed is quite critical, but still without penalizing much the results accuracy :wink:

Thanks again for your time and support!

P.S.: just for asking, do you maybe have a rough idea about when (in which Orekit release) the open bugs on DSST propagator will be fixed? :slight_smile:

5

Hi again @nick,

Nope… :frowning: I’d like to fix it because it annoys me too in my applications, but it’s not easy and I couldn’t find the time yet. Contribution welcome! :wink:

As a general rule, you’ll find that DSST is much much faster in MEAN elements. So what I tend to do is propagate in mean elements and then compute the osculating elements only when I need them. That being said, if you need a lot of osculating elements, like for example an ephemeris with a regular time step, then I suggest you use the numerical propagator, it is more accurate and still quite performant.

Oh, and the SRP in DSST is notoriously slow too, so if it’s not one of the dominant perturbations in your dynamics, I suggest you don’t use as a first-line test.

Cheers,
Maxime

1 Like