Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

New effects to Reboundx - Type I migration and inner disk edge #69

Merged
merged 31 commits into from Feb 23, 2022

Conversation

Kaltrinak
Copy link
Contributor

@Kaltrinak Kaltrinak commented Nov 18, 2021

The reason for this pull request is to contribute a code that implements Type I migration and an inner disk edge to a planetary system.

I am a master student that has been working with mean motion resonance in a two planet system migrating via Type I migration. My master thesis looks at which parameters and how they impact mean motion resonance capture. For this work I needed to write a code for Type I migration as a new effect using the instructions in the Reboundx API and guidance from my supervisor Antoine Petit. In order to stop the planets from falling into the star and since Type I migration is believed to be sequential, an inner disk edge was also necessary. This was added to the code file of Type I migration (to avoid having to import two separate code files when working with Type I migration specifically) and as a separate force effect to be used with other migration types or just constant migration.

The inner disk edge works as a simple function that reverses the migration as the planet reaches close to the inner disk edge, see Pichierri et al 2018. The equations used for Type I migration are based on Cresswell & Nelson 2008, see more on the documentation of the force effects.

This is my first time coding an effect to Reboundx and doing a pull request, tips on improvements are welcomed

@coveralls
Copy link

Coverage Status

Coverage decreased (-12.1%) to 72.727% when pulling b2c3303 on Kaltrinak:master into 063c676 on dtamayo:master.

@dtamayo
Copy link
Owner

dtamayo commented Nov 19, 2021

This is fantastic, thank you! I will go through this more carefully on Monday and we can get this merged!

@dtamayo
Copy link
Owner

dtamayo commented Nov 22, 2021

What a great and careful pull request, thank you!

In terms of making things modular, it would be nice if this effect only updated the tau parameters, and the modify_orbits_forces added the forces, but we don't have things set up so one can ensure that one effect always gets carried out before another (it just depends on the order that the user put them in), so I like your decision to copy that part and make it a self-enclosed effect. It also makes it simpler to incorporate new effects into REBOUNDx and for people to know which paper to cite...

That said, as the number of effects we add multiplies, it does become harder for a user to find which one they want to use...The inner edge effect seems small enough that maybe we could just add it into modify_orbits_forces directly? If hedge or edge is not NULL, it could call rebx_calculating_planet_trap?

What would you think about adding the typeI effect separately, but add the inner edge functionality to the existing effect? I would be very happy to put you as the authors of modify_orbits_forces and put yours as the implementation paper (there's already an overall REBOUNDx paper anyway!) We could then add the inner edge examples as sections of the modify_orbits_forces examples. I'm more than happy to help as needed.

I noticed that in the inner edge Jupyter notebook, the integration doesn't go long enough to reach the inner edge, so we might want to extend the integration timescale so the effect kicks in like in the Type I example.

I think it would also make sense to move the new effect in effects.rst up from the bottom of the file to just above GR, so it shows up in the automated documentation under the Orbit Modifications heading, along with modify_orbits_forces etc.

@Kaltrinak
Copy link
Contributor Author

Thank you for such a positive respond! I am glad the pull request turned out well.

I did not have any prior experience with reboundx, so when discussing with my supervisor how to code the effects, he suggested to copy the code from another file and just change it. It was much easier for sure!

We both like the idea about the inner disk edge you proposed. Just, it might be more "clean" so to say to put the planet trap function as a separate file than import it into modify_orbits_forces and Antoine also suggested it would be good to add the inner disk edge to the modify_orbits_direct too, then the function is always imported ready to use. The function can than be called from the module when dedge and hedge are not null.
Otherwise the function rebx_calculating_planet_trap can be inserted in both files modify_orbits_forces and modify_orbits_direct as I have done in the TypeImig file. In this case the same function will exist in three files, it might be easier for users, but a bit repetitive. What do you think?

As for the implementation paper, we have yet to write it, but yes it would be very nice to be one of the authors and include our paper as one of the implementation papers for sure, but I leave the documentation and that decision to you.

Sure we can add the example as you suggest and fix the integration time for the inner disk edge example, I noticed that I must also fix the initial position for the planet in the numerical simulation to be at 1AU (I had forgotten to change that in the end) so that it is the same as the analytical case otherwise it looks as though the two methods do not match, which they should when the setup is the same.

Sure that does make sense!

@dtamayo
Copy link
Owner

dtamayo commented Nov 23, 2021

That makes sense. Thinking about it more, I wonder whether it would be worthwhile to put in a little bit more work to make this more flexible, given that someone might later want to implement Type II migration or some other disk model, in which case the code duplication could start to get out of hand.

The Type I code is essentially code to update the tau_a, tau_e and tau_inc every timestep. So what about adding a parameter to modify_orbits_forces that's a function pointer, which we could call update_taus. Then you could do

'''
mof = reba.load_force("modify_orbits_forces")
mof.params['update_taus'] = "Type I"
'''

and behind the scenes that sets a function pointer to your typeI function, except it only updates all the taus, and doesn't copy-paste the force calculation. That way if some paper comes along later and wants to advocate for a different prescription, they can

mof.params['update_taus'] = "super fancy Type I" 

and then it's really easy to take your code structure and change the prescription, without duplicating any code.

What do you think? Any opinions @hannorein ? If that sounds reasonable, I'm happy to take your code and write a quick version

@hannorein
Copy link
Collaborator

I can see arguments for both ways. The function pointer is definitely neat. I really like it. But I think a newcomer would have an easier time simply modifying (or copying) an existing effect rather than trying to figure out how the function pointers work.

@acpetit
Copy link
Contributor

acpetit commented Nov 23, 2021

I really like the idea of the update_taus function! In particular, I think that our inner edge prescription is very rough here so it may be good to allow for more physically motivated forms.
What about implementing all the customs effects as update_taus while providing templates within reboundx for the stuff coded by Kaltrina?
I guess one issue becomes how to handle the new parameters (like hedge/dedge for the inner disk example) without having to add them to an evergrowing list.

@dtamayo
Copy link
Owner

dtamayo commented Nov 24, 2021

Decided I'd go for it, and was excited to code something again, but realized that an additional downside is that function pointers would have to be reset when you reload a REBOUNDx binary.
So an advantage of the current implementation is that you can save a binary, reload it, and everything works fine.

I do agree with Hanno that the more fancy things we put in, the less inclined someone new is to try to figure out what everything is doing to write their own effect. It also becomes a bit simpler for each paper that contributes its own effect to be a standalone thing, so it's clear who to cite.

So maybe we could back to the original plan. I like your idea of having a single inner edge file. I'll be around very sporadically the next couple days with the US Thanksgiving holiday, but I'm happy to help as needed!

@acpetit
Copy link
Contributor

acpetit commented Nov 24, 2021

Ok, I think we can handle modifying modify_orbit_forces/direct as well as type_I on our own by adding an update_tau_a function triggered if a flag is provided (as it was suggested initially). I can then add another file for the planet trap so that it doesn't stay into the modify_orbit files. However, I am not sure Kaltrina nor me are familiar enough with the code to provide the framework with the pointer to an update_tau function.

If so do you have any views on the form of the prototype for the function update_tau?
I was thinking about something of the form:
rebx_update_tau(struct reb_simulation* const sim, struct rebx_force* const force, struct reb_particle* p, struct reb_particle* source, double* invtau_a, double* tau_e, double* tau_inc)

@dtamayo
Copy link
Owner

dtamayo commented Nov 25, 2021 via email

@acpetit
Copy link
Contributor

acpetit commented Dec 6, 2021

Hey, I just realized that I haven't give any news in two weeks... I am very busy until New Year but I will take care of this in January!

@dtamayo
Copy link
Owner

dtamayo commented Dec 6, 2021

Sounds great! Thanks Antoine and Kaltrina

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

5 participants