MESA

Modules for Experiments
in Stellar Astrophysics

Latest News

This web documentation corresponds to MESA release r7624.

This page documents the MESA options that are part of the binary_controls namelist. It is autogenerated from the file $MESA_DIR/binary/defaults/binary_controls.defaults.

The documented defaults files available for this version are:

Boxes like

option_name = 'default'

show the default value of each option. To override the default values, add an entry to the binary_controls namelist in your inlist.

specifications for starting model

m1

Initial mass of star 1 in Msun units. Not used when loading a saved model. same caveats as initial_mass in star/defaults/controls.defaults apply

m1 = 1.0d0

m2

Initial mass of star 2 in Msun units. Not used when loading a saved model. same caveats as initial_mass in star/defaults/controls.defaults apply

m2 = 0.8d0

initial_period_in_days

Initial orbital period in days.

initial_period_in_days = 0.5d0

initial_separation_in_Rsuns

Initial separation measured in Rsuns. Only used when initial_period_in_days < 0

initial_separation_in_Rsuns = 100

initial_eccentricity

Initial eccentricity of the system

initial_eccentricity = 0.0

controls for output

history_name

Name of file for binary output

history_name = 'binary_history.data'

append_to_star_history

If true, then the columns from the binary_history are also included in each of the stars history files. NOTE: if .false., then pgstar cannot access binary data

append_to_star_history = .true.

log_directory

Directory for binary output

log_directory = '.'

history_dbl_format

history_int_format

history_txt_format

Format for double, int and text in binary output

history_dbl_format = '(1pes32.16e3, 1x)'
history_int_format = '(i32, 1x)'
history_txt_format = '(a32, 1x)'

photostep

photo_digits

These overwrite the numbers set for photostep and photo_digits for each star, so that profiles are outputted simultaneously

photostep = 50
photo_digits = 3

timestep controls

The terminal output during evolution includes a short string for the 'dt_limit'. This is to give you some indication of what is limiting the time steps. Here's a dictionary mapping those terminal strings to the corresponding control parameters. These only include limits from binary, to see the limits from star refer to star/default/controls.defaults

      terminal output       related parameter
     'b_companion'          timestep limited by companion
     'b_RL'                 fr
     'b_jorb'               fj
     'b_envelope'           fm
     'b_separation'         fa
     'b_eccentricity'       fe
     'b_deltam'             fdm

fm

fa

fr

fj

fe

Timestep controls based on relative changes. After each step an upper limit is set on the timestep based on changes on different quantities. If the quantity is X and the change in one timestep dX, then this limit is given by

dt_next_max = dt * fX*abs(X / dX)

each of these controls deals with the following:

fm = 0.01d0
fa = 0.01d0
fr = 0.10d0
fj = 0.001d0
fe = 0.01d0

fm_limit

fr_limit

fe_limit

Limits to timestep controls give by fm, fr and fe. As these three quantities evolve naturally to zero, following strictly the timestep limit given by fX would reduce timesteps infinetely. These fX_limit avoid this problem by computing the limit to the timestep as

dt_next_max = dt * fX*abs(max(abs(X),fX_limit) / dX)

If any of these fX_limit is smaller than zero it is ignored.

fm_limit = 1d-3
fr_limit = 1d-2
fe_limit = 1d-1

fr_dt_limit

Minimum timestep limit allowed for the fr control in years.

fr_dt_limit = 10d0

fdm

Limits the timestep based on the expected fractional change of mass of the donor,

dt_next_max = fdm * s% mstar / abs(s% star_mdot)
fdm = 0.005d0

dt_softening_factor

Weight factor to average max_timestep with old one (in log space) as in

dt_next_max = 10**(dt_softening_factor*log10(dt_next_max_old) + &
    (1-dt_softening_factor)*log10(dt_next_max))

where dt_next_max_old is the limit used in the previous step. This is meant to avoid large changes in dt. Values must be < 1 and >= 0.

dt_softening_factor = 0.5

varcontrol_{stage}

Allows binary to set varcontrol_target for each star depending on the stage of evolution. Ignored if < 0. Each one controls the following stages,

varcontrol_case_a = -1d0
varcontrol_case_b = -1d0
varcontrol_ms = -1d0
varcontrol_post_ms = -1d0

when to stop

accretor_overflow_terminate

terminate evolution if (r-rl)/rl is bigger than this for accretor

accretor_overflow_terminate = 0.0d0

terminate_if_initial_overflow

terminate evolution if first model of run is overflowing

terminate_if_initial_overflow = .true.

mass transfer controls

mass_transfer_*

Transfer efficiency controls. alpha, beta, delta and gamma parameters as described in Tauris & van den Heuvel 2006 section 16.4.1, transfer efficiency is given by 1-alpha-beta-delta.

These only affect mass that is lost from the donor due to mass transfer, winds from each star will carry away angular momentum from the vicinity of each even when transfer efficiency is unity. Each of these represent the following:

mass_transfer_alpha = 0.0d0
mass_transfer_beta = 0.0d0
mass_transfer_delta = 0.0d0
mass_transfer_gamma = 0.0d0

limit_retention_by_mdot_edd

Limit accretion into the donor using mdot_edd

limit_retention_by_mdot_edd = .true.

use_es_opacity_for_mdot_edd

If .true., then the opacity for mdot_edd is computed as 0.2*(1+X) If .false., the opacity of the outermost cell of the donor is used

use_es_opacity_for_mdot_edd = .false.

use_this_for_mdot_edd

Fixed mdot_edd in Msun/yr, ignored if negative

use_this_for_mdot_edd = -1

mdot_scheme

How to compute mass transfer. Options are:

mdot_scheme = 'Ritter'

explicit mass transfer computation.

MESA can compute mass transfer rates either explicitly (at the beggining of the step) or implicitly (iterating the solution until the mass transfer rate matches the value computed at the end of the step). The explicit method is used if max_tries_to_achieve <= 0.

cur_mdot_frac

Average the explicit mass transfer rate computed with the old in order to smooth large changes.

mass_transfer = mass_transfer_old * cur_mdot_frac + (1-cur_mdot_frac) * mass_transfer
cur_mdot_frac = 0.5d0

max_explicit_abs_mdot

Limit the explicit mass transfer rate to max_explicit_abs_mdot, in Msun/secyer

max_explicit_abs_mdot = 1d-7

implicit mass transfer computation.

max_tries_to_achieve

The implicit method will modify the mass transfer rate and redo the step until it either finds a solution, or the number of tries goes above max_tries_to_achieve. if max_tries_to_achieve <= 0 the explicit method is used.

max_tries_to_achieve = 0

implicit_scheme_tolerance

Tolerance for which a solution is considered valid. For the Ritter and Kolb schemes if we call mdot the mass transfer rate used for the step, and mdot_end the one computed at the end of it, a solution is valid if

|(mdot-mdot_end)/mdot_end| < b% implicit_scheme_tolerance

For the roche_lobe scheme, a solution will be considered valid if

-implicit_scheme_tolerance < (r-rl)/rl < 0

When using the roche_lobe scheme smaller values of order 1d-3 or smaller are recommended.

implicit_scheme_tolerance = 5d-2

initial_change_factor

change_factor_fraction

implicit_lambda

The implicit scheme works by adjusting the mass transfer rate from the previous step until it finds a solution. If the mass transfer needs to increase/reduce after a try, then it is multiplied/divided by change_factor. initial_change_factor provides the initial value for this parameter, however, since at certain points the mass transfer rate will increase steeply and at others remain mostly constant from step to step, MESA adjusts the value of the change factor to make it easier to find solutions. Whenever the mass transfer rate changes from the previous value, MESA will modify the change_factor according to:

if(mass_transfer_rate < mass_transfer_prev) then
   change_factor = change_factor*(1.0-implicit_lambda) &
      + implicit_lambda*(1+change_factor_fraction*(mass_transfer_rate/mass_transfer_prev-1))
else
   change_factor = change_factor*(1.0-implicit_lambda) &
      + implicit_lambda*(1+change_factor_fraction*(mass_transfer_prev/mass_transfer_rate-1))
   change_factor = change_factor*(1.0-implicit_lambda) &
      + implicit_lambda*(1+change_factor_fraction*(mass_transfer_rate/mass_transfer_prev-1))
end if

Choosing implicit_lambda = 0 will keep the change factor constant.

initial_change_factor = 1.5d0
change_factor_fraction = 0.9d0
implicit_lambda = 0.25d0

max_change_factor

min_change_factor

Maximum and minimum values for the change_factor

max_change_factor = 1.5d0
min_change_factor = 1.05d0

starting_mdot

When using the roche_lobe scheme, if the donor overflows for the first time use starting_mdot (in Msun/secyer) as an initial guess for the mass transfer rate.

starting_mdot = 1d-12

roche_min_mdot

When using the roche_lobe scheme, if mass transfer rate is below roche_min_mdot (in Msun/secyer) and the donor is not overflowing its roche lobe, assume detachment and stop mass transfer.

roche_min_mdot = 1d-16

min_mdot_for_implicit

For any choice except for the roche_lobe scheme mass transfer will be computed explicitly until the explicit computation of mdot is > min_mdot_for_implicit (in Msun/secyer), even if max_tries_to_achieve > 0. This is to avoid spending many iterations when the stars are detached and the explicit calculation gives very low values of mdot.

min_mdot_for_implicit = 1d-16

max_implicit_abs_mdot

Limit the implicit mass transfer rate to max_implicit_abs_mdot, in Msun/secyer

max_implicit_abs_mdot = 1d99

Tidal wind enhancement

do_enhance_wind_*

Use the Tout & Eggleton mechanism to tidally enhance the wind mass loss from one or both components according to:

Mdot_w = Mdot_w * ( 1 + B_wind * min( (R/RL)^6, 0.5^6 ) )

Tout & Eggleton 1988,MNRAS,231,823 (eq. 2)

"_1" refers to first star, "_2" to the second one.

do_enhance_wind_1 = .false.
do_enhance_wind_2 = .false.

tout_B_wind_*

The B_wind parameter from the previous equation. Default value is taken from Tout & Eggleton 1988,MNRAS,231,823

"_1" refers to first star, "_2" to the second one.

tout_B_wind_1 = 1d4
tout_B_wind_2 = 1d4

Wind mass accretion

do_wind_mass_transfer_*

transfer part of the mass lost due to stellar winds from the mass losing component to its companion. Using the Bondy-Hoyle mechanism. "_1" refers to first star, "_2" to the second one.

do_wind_mass_transfer_1 = .false.
do_wind_mass_transfer_2 = .false.

wind_BH_alpha_*

Bondy-Hoyle accretion parameter for each star. The default is 3/2 taken from Hurley et al. 2002, MNRAS, 329, 897, in agreement with Boffin & Jorissen 1988, A&A, 205, 155 "_1" refers to first star, "_2" to the second one.

wind_BH_alpha_1 = 1.5d0
wind_BH_alpha_2 = 1.5d0

max_wind_transfer_fraction_*

Upper limit on the wind transfer fraction for star * "_1" refers to first star, "_2" to the second one.

max_wind_transfer_fraction_1 = 0.5d0
max_wind_transfer_fraction_2 = 0.5d0

orbital jdot controls

do_jdot_gr

Include gravitational wave radiation in jdot

do_jdot_gr = .true.

do_jdot_ml

Include loss of angular momentum via mass loss. The parameters mass_transfer_* determine the fractions of mass lost from the vincinity of the donor, the accretor, or a circumbinary coplanar toroid.

do_jdot_ml = .true.

do_jdot_ls

Fix jdot such that the total angular momentum of the system is conserved, except for loses due to other jdot mechanisms, or angular momentum loss from winds. This is meant to take care of L-S coupling due to tides.

do_jdot_ls = .true.

do_jdot_missing_wind

Usually MESA computes stellar AM loss due to winds by taking the angular momentum from the removed layers of the star. However, when mass transfer is included, wind mass loss and mass accretion are added up, and only the remainder, if corresponding to net mass loss, contributes to stellar AM loss. jdot_missing_wind compensates for this, by removing from the orbit an amount of angular momentum equal to the mass lost that does not contribute to stellar AM loss, times the specific angular momentum at the surface.

do_jdot_missing_wind = .false.

do_jdot_mb

Include magnetic braking as in Rappaport, Verbunt, and Joss. apj, 275, 713-731. 1983.

do_jdot_mb = .true.

include_accretor_mb

If true, the contribution to jdot from magnetic braking of the accretor is also taken into account.

include_accretor_mb = .false.

magnetic_braking_gamma

gamma exponent for magnetic braking.

magnetic_braking_gamma = 3.0d0

keep_mb_on

If true keep magnetic braking even when radiative core goes away.

keep_mb_on = .false.

jdot_multiplier

Multiply total jdot by this factor. NOTE: jdot_ls is not affected by this.

jdot_multiplier = 1d0

rotation and sync controls

do_j_accretion

If true, compute accretion of angular momentum following A.3.3 of de Mink et al. 2013, ApJ, 764, 166. Otherwise, incoming material is assumed to have the specific angular momentum of the surface of the accretor.

do_j_accretion = .false.

do_tidal_sync

If true, apply tidal torque to the star

do_tidal_sync = .false.

sync_type_*

Timescale for orbital synchronisation. "_1" refers to first star, "_2" to the second one. Options are:

sync_type_1 = 'Hut_conv'
sync_type_2 = 'Hut_conv'

sync_mode_*

Where angular momentum is deposited for synchronization. "_1" refers to first star, "_2" to the second one. Options are:

sync_mode_1 = 'Uniform'
sync_mode_2 = 'Uniform'

Ftid_*

Tidal strength factor. Synchronisation and circularisation timescales are divided by this. "_1" refers to first star, "_2" to the second one.

Ftid_1 = 1d0
Ftid_2 = 1d0

do_initial_orbit_sync_*

Relax rotation of star to orbital period at the beggining of evolution. "_1" refers to first star, "_2" to the second one.

do_initial_orbit_sync_1 = .false.
do_initial_orbit_sync_2 = .false.

tidal_reduction

tidal_reduction accounts for the reduction in the effectiveness of convective damping of the equilibrium tide when the tidal forcing period is less than the convective turnover period of the largest eddies. It corresponds to the exponent in eq. (32) of Hurley et al. 2002, MNRAS, 329, 897

tidal_reduction = 1 follows Zahn(1966, 1989), while tidal_reduction = 2 follows Goldreich & Nicholson (1977).

tidal_reduction = 2.0d0

eccentricity controls

do_tidal_circ

If true, apply tidal circularisation

do_tidal_circ = .false.

circ_type_*

Mechanism for circularisation. Options are: "_1" refers to first star, "_2" to the second one.

circ_type_1 = 'Hut_conv'
circ_type_2 = 'Hut_conv'

use_eccentricity_enhancement

Flag to turn on Soker eccentricity enhancement

use_eccentricity_enhancement = .false.

max_abs_edot_tidal

Maximum absolute value for tidal edot (in 1/s). If the computed tidal edot goes above this, then it is fixed at this maximum

max_abs_edot_tidal = 1d-6

max_abs_edot_enhance

Maximum absolute value for eccentricity enhancement (in 1/s). If the computed edot goes above this, then it is fixed at this maximum

max_abs_edot_enhance = 1d-6

min_eccentricity

If after a step eccentricity < min_eccentricity, then fix it at this value

min_eccentricity = 0.0

max_eccentricity

If after a step eccentricity > max_eccentricity, then fix it at this value

max_eccentricity = 0.99

irradiation controls

accretion_powered_irradiation

Flag to turn on irradiation of the donor due to accretion onto a compact object.

accretion_powered_irradiation = .false.

accretor_radius_for_irrad

Radius in cm for accreting object, e.g., 10km for ns.

accretor_radius_for_irrad = 1d6

col_depth_for_eps_extra

Energy from irradiation will be deposited in the outer 4*Pi*R^2*col_depth_for_eps_extra grams of the star.

col_depth_for_eps_extra = -1

irrad_flux_at_std_distance

std_distance_for_irradiation

If irrad_flux_at_std_distance > 0 then irradiation flux is computed as

s% irradiation_flux = b% irrad_flux_at_std_distance * &
    (b% std_distance_for_irradiation/b% separation)**2
irrad_flux_at_std_distance = -1
std_distance_for_irradiation = -1

max_F_irr

Limit irradiation by this amount.

max_F_irr = 5d12

miscellaneous controls

keep_donor_fixed

keep star 1 as donor, even if accretor is closer to filling roche lobe

keep_donor_fixed = .true.

mdot_limit_donor_switch

Do not change donor if mass transfer is larger than this (given in Msun/secyer). Avoids erratic changes when both stars are filling their roche loches.

mdot_limit_donor_switch = 1d-20

use_other_{hook}

Logicals to deploy the use_other routines.

use_other_rlo_mdot = .false.
use_other_tsync = .false.
use_other_mdot_edd = .false.
use_other_accreted_material_j = .false.
use_other_jdot_gr = .false.
use_other_jdot_ml = .false.
use_other_jdot_ls = .false.
use_other_jdot_missing_wind = .false.
use_other_jdot_mb = .false.
use_other_extra_jdot = .false.
use_other_binary_wind_transfer = .false.