hannorein / rebound

💫 An open-source multi-purpose N-body code.
https://rebound.readthedocs.io/
GNU General Public License v3.0
861 stars 221 forks source link

integrating simulation and extract ephemerides at a given JD #564

Closed ifermion closed 3 years ago

ifermion commented 3 years ago

Hi all,

First of all, let me start off by thanking you all for such a amazing tool that rebound is, while I'm still discovering it.

Let's say you import the ephemerides from JPL using sim.add() and using the argument date to specify the Julian Date (JD) at which you want to instantiate the simulation (sim.t = 0). You keep integrating backwards (with something like sim.dt = -0.001) for let's say 100 years (preferably using 100*2pi, where 2pi is a year as I've seen in the examples). Then, say you have a list of specific JDs/ epochs that correspond to this period of 100 years, for which you want to extract ephermerides for certain objects in your simulation (like distance between two planets for example) as the integrations runs.

So I'm wondering how one would go about this task. Any ideas ? There is a way to translate sim.t to JD ?

Please let me know if something is not clear, I'll try to explain better.

Thanks in advance !

Rmelikyan commented 3 years ago

Hi,

From what I've seen, tracking time is not a native feature to rebound and such I'm sure each of us have our own solutions for corresponding rebounds internal clock with the physical systems we may be representing.

One thing you can do (what I do at least) is use pythons datetime object to keep track of the initial time of a simulation. You can then use the timedelta object to modify the initial time by some number of seconds or days. This might look like:

t_initial = datetime.datetime(2021, 8, 25, 12, 10, 0) sim.integrate(1000) # Here I'm assuming you are using SI units so that 1000 coresponds to 1000 seconds. If that is not the case then you may want to do some conversions in the next step new_time = t_initial + datetime.timedelta(d=0, s=1000) # here I make a new datetime object that is 1000 seconds later than the initial time and therefore matches the simulation.

@hannorein If I find time I'll try to put together a small example notebook that looks a little bit deeper into time tracking in rebound and what options you may have available...

-Robert

hannorein commented 3 years ago

Hi,

Sorry for the slow response.

It's pretty straightforward. All you need to know is that one Julian year is defined as 365.25 days and that one REBOUND time unit is by default equal to one year/(2*pi) ... (because G=1 in code units).

Say you import the ephemeris at some Julian day JD_0 and set the simulation time (sim.t) to 0 at the beginning. Then, you can integrate for some time T in REBOUND. In Julian days, that will correspond to JD_0 + T/(2.pi)365.25. That'll be accurate to about 1e-5 after 1 year. If you want it more accurate you need to take into account the difference between a Julian year and an actual year.

(Robert posted some helpful snippets if you want the Julian date in days, months, years rather than the standard decimal format)

I hope that helps.

Hanno

danielk333 commented 3 years ago

This is just in addition to what @hannorein suggested above. For converting between different times I usually use astropy.time. Especially since its allows seamless conversion between TAI (International Atomic Time) and UTC which is often used by instruments (radars ect) and TDB (Barycentric Dynamical Time) used by JPL Kernels and a few other sources.

ifermion commented 3 years ago

Hi @Rmelikyan, @hannorein, @danielk333

Thank you very much for your replies. After some testing, I think the approch of using the python datetime object and swtiching to SI units described by @Rmelikyan appears the most straighforward to me. I'm also using astropy.time Time objects to swiftly convert between UTC and JD.

Best regards !