dtamayo / reboundx

A library for adding additional forces to the REBOUND N-body integration package
GNU General Public License v3.0
80 stars 60 forks source link

Add Tidal Drag Effect #44

Closed sabaronett closed 4 years ago

sabaronett commented 4 years ago

Aside from:

I've completed the "Adding A New Effect" guide for the first tides_drag implementation and am submitting the first pull request for review.

coveralls commented 4 years ago

Coverage Status

Coverage increased (+3.3%) to 73.384% when pulling a00c4467774bc85ac2b4fd71086bff11a8f74001 on sabaronett:tides-drag into 39320c94ebc705a801d3dd50596fbe86e04b2d96 on dtamayo:master.

codecov-io commented 4 years ago

Codecov Report

:exclamation: No coverage uploaded for pull request base (master@27dab18). Click here to learn what that means. The diff coverage is 50%.

Impacted file tree graph

@@            Coverage Diff            @@
##             master      #44   +/-   ##
=========================================
  Coverage          ?   72.85%           
=========================================
  Files             ?       38           
  Lines             ?     3408           
  Branches          ?        0           
=========================================
  Hits              ?     2483           
  Misses            ?      925           
  Partials          ?        0
Impacted Files Coverage Δ
src/track_min_distance.c 0% <ø> (ø)
src/central_force.c 82.53% <ø> (ø)
src/radiation_forces.c 80% <ø> (ø)
src/gravitational_harmonics.c 98.43% <ø> (ø)
src/tides_drag.c 0% <ø> (ø)
src/core.c 82.5% <ø> (ø)
src/modify_orbits_forces.c 97.43% <ø> (ø)
reboundx/extras.py 91.01% <100%> (ø)
reboundx/params.py 95.23% <40%> (ø)

Continue to review full report at Codecov.

Legend - Click here to learn more Δ = absolute <relative> (impact), ø = not affected, ? = missing data Powered by Codecov. Last update 27dab18...a00c446. Read the comment docs.

sabaronett commented 4 years ago

@dtamayo, I just made some minor modifications (mostly better comments and freeing memory at the end) to example/problem.c in 761387fafeade6163b02da5489eb185cf528054f on my repo's tides_drag branch. I'm not quite sure how that exactly plays out with this in-progress pull-request...

sabaronett commented 4 years ago

Never mind! I just checked the details under the "Files changed" tab at the top of this request and verified my latest commits are there. Sorry to bother!

dtamayo commented 4 years ago

So I spent some time going through this, going through the literature, and doing some soul searching! The issue I see with tides (and with a library of effects more generally) is that everyone has a different favorite prescription. That becomes very unwieldy if everyone adds their own as a separate effect. That's confusing with different registered param names etc., and even just for a user looking through the effects. I think the best way forward is to bundle effects by their functional form, and parametrize the overall factors. So a tidal effect might have a Q or a time lag. Then when someone comes along with a new favorite parametrization, instead of writing a new effect, they write a new convenience function that calculates the appropriate Q or time lag from their favorite set of parameters. For example, here, we could take Zahn's formalism with convective timescales etc. and write a convenience function to calculate the corresponding time lag.

Going through the equations I also realized that we can make your effect much more general (have it work with eccentric orbits and provide eccentricity damping) with a couple small changes. How would you feel about combining this with tides_precession (which gives the radial component)? That would give us a general implementation of the constant time lag tidal model. I would be happy to help implement the required changes, and we would still make you the lead author of the combined effect.

sabaronett commented 4 years ago

Thanks for spending the time to go through everything in detail. I completely agree with your outlook to avoid what could end up becoming a hodgepodge of user-submitted effects, with some being only slight variations of others. And as you point out, this effect would indeed serve better if generalized, rather than have a limited application to just circular orbits. Combining it with tides_precession is perfectly fine. Please just let me know how I can help with that effort in any way.

Given this new development, should I hold off then on looking into and devising an analytical solution test? At least until its more final implementation within tides_precession? Let me know, as I can instead work on the cubic spline interpolation when I have time again to work on this project.

dtamayo commented 4 years ago

Awesome. Don't worry about the analytical test. I will implement and test. I'll update you here and maybe you could add a convenience function to calculate the relevant tidal parameters from the constants you're using from Zahn?

sabaronett commented 4 years ago

Sure thing. I can have them ready by next week. Let me know if you need them done sooner.

sabaronett commented 4 years ago

I'm not sure if I got exactly what you were looking for, but I added some convenience functions to calculate the tidal parameters in a00c4467774bc85ac2b4fd71086bff11a8f74001: https://github.com/sabaronett/reboundx/blob/a00c4467774bc85ac2b4fd71086bff11a8f74001/src/tides_drag.c#L129-L147

Please let me know if I'm on the right track, and I apologize if I didn't quite understand and missed the mark.

dtamayo commented 4 years ago

Hi Stanley,

I've implemented and tested the constant time lag model now on the tides_drag branch of my reboundx repository. Could you check out the TidesConstantTimeLag ipytthon example and check that it does what you think it should? I've added several unit tests that seem to show it's consistent with what you were getting with the previous implementations, and it works in the more. general case now.

What's left to do is updating the C example and the documentation. If you see any errors/typos along the way please feel free to fix them! Once you're happy we can merge this effect

sabaronett commented 4 years ago

Hi Dan,

Just finished going through the Jupyter notebook, and it looks pretty good. Other than these output warnings:

/Users/dtamayo/miniconda3/envs/p3/lib/python3.7/site-packages/matplotlib/font_manager.py:1331: UserWarning: findfont: Font family ['serif'] not found. Falling back to DejaVu Sans
  (prop.get_family(), self.defaultFamily[fontext]))

I only found one minor typo (missing apostrophe):

We additionally assume that bodies['] spins are aligned with the reference z axis.

I would fix it myself, but I'm not sure how since it's on your repo's branch.

Let me know about the C example or documentation, so I can review/proofread when ready, or if you'd like me to help with either.

One question I had about the notebook was cell [6]:

apred = ps[0].r*((ps[1].a/ps[0].r)**8 - 48.*ps[0].params["tctl_k1"]*q*(1+q)*times/T)**(1./8.)

I see where you express T from Hut Eq. 12. But I can't seem to make the connection between apred and Eq. 9. Is there something else I'm missing, or a different equation perhaps?

I appreciate what you've done in generalizing this to eccentric orbits. Maybe in future releases, we can further generalize the effect to include inclined orbits based on the rest of Hut's paper!

sabaronett commented 4 years ago

Actually, I suppose I could fork the notebook to a new branch on my repo, make the corrections, then submit a pull request. Let me know if you'd like that.

dtamayo commented 4 years ago

Thanks Stanley, will fix. Yeah I agree it's opaque...I solved Hut's diff eq for da/dt assuming e=0 two get apred.

And the current implementation also works for inclined orbits. The major limitation is not self-consistently evolving the spins, but that's a fairly major addition that I can't work on now. Maybe having a starting point will motivate someone that needs it to include it!

Don't worry. I'll go through it and merge, we can improve clarity and fix any typos after the fact. I wanted to make sure all tests passed

dtamayo commented 4 years ago

I've now merged the other branch including all your commits. Thanks so much for all the hard work, and for all the great fixes to the documentation! Check it out. Closing this now.