Skip to content

Commit

Permalink
Merge b2c3303 into 063c676
Browse files Browse the repository at this point in the history
  • Loading branch information
@Kaltrinak
Kaltrinak committed Nov 18, 2021
2 parents 063c676 + b2c3303 commit 8441281
Show file tree
Hide file tree
Showing 14 changed files with 1,000 additions and 47 deletions.
7 changes: 7 additions & 0 deletions doc/effect_headers.txt
View file
Expand Up @@ -57,3 +57,10 @@ $$$$$$$$$$$$
Miscellaneous Utilities
^^^^^^^^^^^^^^^^^^^^^^^

$$$$$$$$$$$$$$$$$
Inner disk edge
^^^^^^^^^^^^^^^^^

$$$$$$$$$$$$$$$$$$
Type I migration
^^^^^^^^^^^^^^^^^^
187 changes: 143 additions & 44 deletions doc/effects.rst
View file
Expand Up @@ -16,6 +16,50 @@ Positive timescales correspond to growth / progression, negative timescales corr
Semimajor axes, eccentricities and inclinations grow / damp exponentially.
Pericenters and nodes progress/regress linearly.

.. _modify_orbits_forces:

modify_orbits_forces
********************

======================= ===============================================
Authors D. Tamayo, H. Rein
Implementation Paper `Kostov et al., 2016 <https://ui.adsabs.harvard.edu/abs/2016ApJ...832..183K/abstract>`_.
Based on `Papaloizou & Larwood 2000 <http://labs.adsabs.harvard.edu/adsabs/abs/2000MNRAS.315..823P/>`_.
C Example :ref:`c_example_modify_orbits`
Python Example `Migration.ipynb <https://github.com/dtamayo/reboundx/blob/master/ipython_examples/Migration.ipynb>`_
`EccAndIncDamping.ipynb <https://github.com/dtamayo/reboundx/blob/master/ipython_examples/EccAndIncDamping.ipynb>`_.
======================= ===============================================

This applies physical forces that orbit-average to give exponential growth/decay of the semimajor axis, eccentricity and inclination.
The eccentricity damping keeps the angular momentum constant (corresponding to `p=1` in modify_orbits_direct), which means that eccentricity damping will induce some semimajor axis evolution.
Additionally, eccentricity/inclination damping will induce pericenter/nodal precession.
Both these effects are physical, and the method is more robust for strongly perturbed systems.

**Effect Parameters**

If coordinates not, defaults to using Jacobi coordinates.

============================ =========== ==================================================================
Field (C type) Required Description
============================ =========== ==================================================================
coordinates (enum) No Type of elements to use for modification (Jacobi, barycentric or particle).
See the examples for usage.
============================ =========== ==================================================================

**Particle Parameters**

One can pick and choose which particles have which parameters set.
For each particle, any unset parameter is ignored.

============================ =========== ==================================================================
Field (C type) Required Description
============================ =========== ==================================================================
tau_a (double) No Semimajor axis exponential growth/damping timescale
tau_e (double) No Eccentricity exponential growth/damping timescale
tau_inc (double) No Inclination axis exponential growth/damping timescale
============================ =========== ==================================================================


.. _exponential_migration:

exponential_migration
Expand Down Expand Up @@ -93,50 +137,6 @@ tau_omega (double) No Period of linear apsidal precession/reg
============================ =========== ==================================================================


.. _modify_orbits_forces:

modify_orbits_forces
********************

======================= ===============================================
Authors D. Tamayo, H. Rein
Implementation Paper `Kostov et al., 2016 <https://ui.adsabs.harvard.edu/abs/2016ApJ...832..183K/abstract>`_.
Based on `Papaloizou & Larwood 2000 <http://labs.adsabs.harvard.edu/adsabs/abs/2000MNRAS.315..823P/>`_.
C Example :ref:`c_example_modify_orbits`
Python Example `Migration.ipynb <https://github.com/dtamayo/reboundx/blob/master/ipython_examples/Migration.ipynb>`_
`EccAndIncDamping.ipynb <https://github.com/dtamayo/reboundx/blob/master/ipython_examples/EccAndIncDamping.ipynb>`_.
======================= ===============================================

This applies physical forces that orbit-average to give exponential growth/decay of the semimajor axis, eccentricity and inclination.
The eccentricity damping keeps the angular momentum constant (corresponding to `p=1` in modify_orbits_direct), which means that eccentricity damping will induce some semimajor axis evolution.
Additionally, eccentricity/inclination damping will induce pericenter/nodal precession.
Both these effects are physical, and the method is more robust for strongly perturbed systems.

**Effect Parameters**

If coordinates not, defaults to using Jacobi coordinates.

============================ =========== ==================================================================
Field (C type) Required Description
============================ =========== ==================================================================
coordinates (enum) No Type of elements to use for modification (Jacobi, barycentric or particle).
See the examples for usage.
============================ =========== ==================================================================

**Particle Parameters**

One can pick and choose which particles have which parameters set.
For each particle, any unset parameter is ignored.

============================ =========== ==================================================================
Field (C type) Required Description
============================ =========== ==================================================================
tau_a (double) No Semimajor axis exponential growth/damping timescale
tau_e (double) No Eccentricity exponential growth/damping timescale
tau_inc (double) No Inclination axis exponential growth/damping timescale
============================ =========== ==================================================================


General Relativity
^^^^^^^^^^^^^^^^^^

Expand Down Expand Up @@ -515,3 +515,102 @@ min_distance_orbit (reb_orbit) No Parameter to store orbital elements
================================ =========== =======================================================


Inner disk edge
^^^^^^^^^^^^^^^^^

Type I migration
^^^^^^^^^^^^^^^^^^.. _type_I_migration:

type_I_migration
****************

======================= ========================================================================================================================================================================================
Authors Kajtazi, Kaltrina and D. Petit, C. Antoine
Implementation Paper `Kajtazi et al. in prep.
Based on `Cresswell & Nelson 2008 <https://ui.adsabs.harvard.edu/abs/2008A%26A...482..677C/abstract>, and Pichierri et al 2018 <https://ui.adsabs.harvard.edu/abs/2018CeMDA.130...54P/abstract>.
C example :ref: `c_examples_type_I_migration`
Python example `TypeIMigration.ipynb <https://github.com/dtamayo/reboundx/blob/master/ipython_examples/TypeIMigration.ipynb>`_.
======================= ========================================================================================================================================================================================

This applies Type I migration, where eccentricity, semi-major axis and inclination are dampened during migration.
The base of the code is the same as the modified orbital forces one written by D. Tamayo, H. Rein.
Moreover, the first part of the code below is the implementation of an inner disc edge, which is decribed and written in the same way in a separate file too,
because it can then be used on its own with another migration precription too, not just in connection with this Type I migration prescription. The inner disc edge is included here directly for
simplicity instead of having to add both separately when using this Type I migration prescription.
Note that this code is not machine independent since power laws were not possible to avoid all together.

**Effect Parameters**

============================ =========== ==================================================================================================================
Field (C type) Required Description
============================ =========== ==================================================================================================================
dedge (double) Yes The position of the inner disk edge in code units
hedge (double) Yes The aspect ratio at the inner disk edge; the disk edge width
sd0 (double) Yes Disk surface density at one code unit from the star; used to find the surface density at any distance from the star
h0 (double) Yes The scale height at one code unit from the star; used to find the aspect ratio at any distance from the star
s (double) Yes Exponent of disk surface density, indicative of the surface density profile of the disk
beta (double) Yes The flaring index; 1 means disk is irradiated by only the stellar flux
============================ =========== ==================================================================================================================

**Particle Parameters**

One can pick and choose which particles have which parameter set.

============================ =========== ===================================================================================
Field (C type) Required Description
============================ =========== ===================================================================================
tau_a (double) No Semimajor axis exponential growth/damping timescale
tau_e (double) No Eccentricity exponential growth/damping timescale
tau_inc (double) No Inclination axis exponential growth/damping timescale
tau_a_red (double) No Planet trap function to stop further migration once the inner disk edge is reached
============================ =========== ===================================================================================


inner_disk_edge
^^^^^^^^^^^^^^^

.. _inner_disk_edge:

inner_disk_edge
***************


$Inner disk edge$ // Effect category

======================= ============================================================================================
Authors Kajtazi, Kaltrina and D. Petit, C. Antoine
Implementation Paper `Kajtazi et al. in prep.
Based on `Pichierri et al 2018 <https://ui.adsabs.harvard.edu/abs/2018CeMDA.130...54P/abstract>`_.
C example :ref: `c_examples_inner_disk_edge`
Python example `InnerDiskEdge.ipynb <https://github.com/dtamayo/reboundx/blob/master/ipython_examples/InnerDiskEdge.ipynb>`_.
======================= ============================================================================================

This applies an inner disk edge that functions as a planet trap. Within its width the planet's migration is reversed
by an opposite and roughly equal magnitude torque. Thus, stopping further migration and trapping the planet within
the width of the trap. The base used here is modified_orbital_forces script written by D. Tamayo, H. Rein.
This implementation should work with any migration/effect not just Type I migration or constant migration.
Other precriptions have not been tested but should work fine, as long as that migration prescription can be given
in terms of the timescales of change in orbital elements and applied through accelerations as done here.

**Effect Parameters**

============================ =========== ===================================================================================
Field (C type) Required Description
============================ =========== ===================================================================================
dedge (double) Yes The position of the inner disk edge in code units
hedge (double) Yes The aspect ratio at the inner disk edge; the disk edge width
============================ =========== ===================================================================================

**Particle Parameters**

One can pick and choose which particles have which parameters set.

============================ =========== ===================================================================================
Field (C type) Required Description
============================ =========== ===================================================================================
tau_a (double) No Semimajor axis exponential growth/damping timescale
tau_e (double) No Eccentricity exponential growth/damping timescale
tau_inc (double) No Inclination axis exponential growth/damping timescale
tau_a_red (double) No Planet trap function to stop further migration once the inner disk edge is reached
============================ =========== ===================================================================================

48 changes: 48 additions & 0 deletions examples/inner_disk_edge/Makefile
View file
@@ -0,0 +1,48 @@
export OPENGL=1

ifndef REB_DIR
ifneq ($(wildcard ../../../rebound/.*),) # Check for REBOUND in default location
REB_DIR=../../../rebound
endif
ifneq ($(wildcard ../../../../rebound/.*),) # Check for REBOUNDx being inside REBOUND directory
REB_DIR=../../../
endif
endif
ifndef REB_DIR # REBOUND is not in default location and REB_DIR is not set
$(error REBOUNDx not in the same directory as REBOUND. To use a custom location, you Must set the REB_DIR environment variable for the path to your rebound directory, e.g., export REB_DIR=/Users/dtamayo/rebound. See reboundx.readthedocs.org)
endif
PROBLEMDIR=$(shell basename `dirname \`pwd\``)"/"$(shell basename `pwd`)

include $(REB_DIR)/src/Makefile.defs

REBX_DIR=../../

all: librebound.so libreboundx.so
@echo ""
@echo "Compiling problem file ..."
$(CC) -I$(REBX_DIR)/src/ -I$(REB_DIR)/src/ -Wl,-rpath,./ $(OPT) $(PREDEF) problem.c -L. -lreboundx -lrebound $(LIB) -o rebound
@echo ""
@echo "Problem file compiled successfully."

librebound.so:
@echo "Compiling shared library librebound.so ..."
$(MAKE) -C $(REB_DIR)/src/
@echo "Creating link for shared library librebound.so ..."
@-rm -f librebound.so
@ln -s $(REB_DIR)/src/librebound.so .

libreboundx.so:
@echo "Compiling shared library libreboundx.so ..."
$(MAKE) -C $(REBX_DIR)/src/
@-rm -f libreboundx.so
@ln -s $(REBX_DIR)/src/libreboundx.so .

clean:
@echo "Cleaning up shared library librebound.so ..."
@-rm -f librebound.so
$(MAKE) -C $(REB_DIR)/src/ clean
@echo "Cleaning up shared library libreboundx.so ..."
@-rm -f libreboundx.so
$(MAKE) -C $(REBX_DIR)/src/ clean
@echo "Cleaning up local directory ..."
@-rm -vf rebound
49 changes: 49 additions & 0 deletions examples/inner_disk_edge/problem.c
View file
@@ -0,0 +1,49 @@
/**
* Inner disk edge.
*
* This example shows how to add an inner disk edge.
*/
#include <stdio.h>
#include <string.h>
#include "rebound.h"
#include "reboundx.h"
#include "core.h"

int main(int argc, char* argv[]){
struct reb_simulation* sim = reb_create_simulation();
/*sim units are ('yr', 'AU', 'Msun')*/
sim->G = 4*M_PI*M_PI;

struct reb_particle star = {0};
star.m = 1.;
reb_add(sim, star);

double m = 0.00001;
double a = 1;
double e = 0;
double inc = 0.;
double Omega = 0.;
double omega = 0.;
double f = 0.;

struct reb_particle p1 = reb_tools_orbit_to_particle(sim->G, star, m, a, e, inc, Omega, omega, f);
reb_add(sim, p1);

sim->dt = 0.002; //The period at inner disk edge divided by 20, for a disk edge location at 0.1 AU
sim->integrator = REB_INTEGRATOR_WHFAST;

struct rebx_extras* rebx = rebx_attach(sim);

struct rebx_force* inner_edge = rebx_load_force(rebx, "inner_disk_edge");
rebx_set_param_double(rebx, &inner_edge->ap, "inner_disk_edge_position", 0.1);
rebx_set_param_double(rebx, &inner_edge->ap, "disk_edge_width", 0.02); //Calculated using a scale height value of 0.03. See Pichierri et al. 2018 for the equation
rebx_add_force(rebx, inner_edge);

double tmax = 1.e4;
rebx_set_param_double(rebx, &sim->particles[1].ap, "tau_a", -tmax);
rebx_set_param_double(rebx, &sim->particles[1].ap, "tau_e", -tmax/100.);

reb_integrate(sim, tmax);
rebx_free(rebx);
reb_free_simulation(sim);
}
48 changes: 48 additions & 0 deletions examples/type_I_migration/Makefile
View file
@@ -0,0 +1,48 @@
export OPENGL=1

ifndef REB_DIR
ifneq ($(wildcard ../../../rebound/.*),) # Check for REBOUND in default location
REB_DIR=../../../rebound
endif
ifneq ($(wildcard ../../../../rebound/.*),) # Check for REBOUNDx being inside REBOUND directory
REB_DIR=../../../
endif
endif
ifndef REB_DIR # REBOUND is not in default location and REB_DIR is not set
$(error REBOUNDx not in the same directory as REBOUND. To use a custom location, you Must set the REB_DIR environment variable for the path to your rebound directory, e.g., export REB_DIR=/Users/dtamayo/rebound. See reboundx.readthedocs.org)
endif
PROBLEMDIR=$(shell basename `dirname \`pwd\``)"/"$(shell basename `pwd`)

include $(REB_DIR)/src/Makefile.defs

REBX_DIR=../../

all: librebound.so libreboundx.so
@echo ""
@echo "Compiling problem file ..."
$(CC) -I$(REBX_DIR)/src/ -I$(REB_DIR)/src/ -Wl,-rpath,./ $(OPT) $(PREDEF) problem.c -L. -lreboundx -lrebound $(LIB) -o rebound
@echo ""
@echo "Problem file compiled successfully."

librebound.so:
@echo "Compiling shared library librebound.so ..."
$(MAKE) -C $(REB_DIR)/src/
@echo "Creating link for shared library librebound.so ..."
@-rm -f librebound.so
@ln -s $(REB_DIR)/src/librebound.so .

libreboundx.so:
@echo "Compiling shared library libreboundx.so ..."
$(MAKE) -C $(REBX_DIR)/src/
@-rm -f libreboundx.so
@ln -s $(REBX_DIR)/src/libreboundx.so .

clean:
@echo "Cleaning up shared library librebound.so ..."
@-rm -f librebound.so
$(MAKE) -C $(REB_DIR)/src/ clean
@echo "Cleaning up shared library libreboundx.so ..."
@-rm -f libreboundx.so
$(MAKE) -C $(REBX_DIR)/src/ clean
@echo "Cleaning up local directory ..."
@-rm -vf rebound

0 comments on commit 8441281

Please sign in to comment.