Problem with the example: full version of yarkovsky_effect: Python 2.7

[PK-ReboundX]

I have a problem with the Yarkovsky effect example: the simplified version (flags: -1,1) works equally and is portable from Python 3.x to Python 2.7. In contrast, the 'full' version of the effect with flag '0' gives different results on Python 2.7 (as if the effect works very poorly or not at all).

All these comments apply to the example included in the Reboundx documentation. Does this mean that I should abandon the version that is currently running on Python 2.7? I work on a remote system and so far everything has worked.

Edit: please note that in all cases I use rebound 3.22.0 and reboundx 3.8.0 (Python 2.7 and Python 3.x: 3.7, 3.8...).

Thanks in advance for any suggestion.

[Nofe4108]

Thanks for reporting this issue. I'll try and figure out if there are any possible solutions as soon as I can. In the meantime, could you please tell me the exact results you're getting from each version of Python.

[PK-ReboundX]

Thanks for response; I will provide more details:

I tried the full version of effect provided in example (link: https://github.com/dtamayo/reboundx/blob/master/ipython_examples/YarkovskyEffect.ipynb)

My results on two computers:

1. [Python 2.7.16, rebound 3.22.0, reboundx 3.8.0]

OUTPUT:

INITIAL ORBITS: <rebound.Orbit instance, a=0.5 e=1.79989376117e-16 inc=0.0 Omega=0.0 omega=0.0 f=0.0> ('CHANGE IN SEMI-MAJOR AXIS:', -6.98876373439461e-09, 'AU \n')

2. [Python 3.8.10, rebound 3.22.0, reboundx 3.8.0]

OUTPUT:

INITIAL ORBITS: <rebound.Orbit instance, a=0.5000000000000001 e=1.799893761170345e-16 inc=0.0 Omega=0.0 omega=0.0 f=0.0> CHANGE IN SEMI-MAJOR AXIS: 1.921943330540632e-05 AU

Please note that the 'simplified version' give the same results of drift. I have a problem with 'full version' only.

Regards, PK

[Nofe4108]

Thanks for sharing this info. I ran the example using Python 2.7 and received the same result you did. Good news is I think I've found the solution to this issue.

In the code, the thermal inertia of the object is given in the following line:

Gamma = (310*u.kg/u.s**(5/2)).to(u.Msun/u.yr**(5/2)) #thermal inertia of object

The issue lies with the exponents. When you perform integer division in Python 2.7, it rounds the result to an integer instead of giving you back a float. In this case, 5/2 was rounded down to 2. Changing this line to something like this:

Gamma = (310*u.kg/u.s**(5.0/2)).to(u.Msun/u.yr**(5.0/2)) #thermal inertia of object

should give you the same answer you get with Python 3.8.

[PK-ReboundX]

This is very good news. Looks like a typical accidental portability error. I changed the same values from integer to floating point and it fixed itself. Thank you for your quick help, I will continue to experiment with this effect.

Regards,

PK

[Nofe4108]

No problem. Please let us know if you run into any other issues.

[dtamayo]

updated that ipython example to also work in python 2.7

[PK-ReboundX]

Please note that one of the fractions in the formula (5/2) was left untouched.

[dtamayo]

thanks for catching that