Planetary Orbital Evolution due to Tides
Orbital evolution of two objects experiencing tides
Evolve::BrokenPowerlawPhaseLagZone Class Reference

A DissipatingZone where the phase lag is described by a broken powerlaw. More...

#include <BrokenPowerlawPhaseLagZone.h>

+ Inheritance diagram for Evolve::BrokenPowerlawPhaseLagZone:
+ Collaboration diagram for Evolve::BrokenPowerlawPhaseLagZone:

Public Member Functions

 BrokenPowerlawPhaseLagZone ()
 Create a non-dissipative zone. Must call setup() if the zone is dissipative. More...
 
void setup (const std::vector< double > &tidal_frequency_breaks, const std::vector< double > &spin_frequency_breaks, const std::vector< double > &tidal_frequency_powers, const std::vector< double > &spin_frequency_powers, double reference_phase_lag, double inertial_mode_enhancement=1.0, double inertial_mode_sharpness=10.0)
 Seup the zone with the given breaks/powers and inertial mode enhancement. Continuous accress all breaks. Must only be called before use. More...
 
virtual void configure (bool initialize, double age, double orbital_frequency, double eccentricity, double orbital_angmom, double spin, double inclination, double periapsis, bool spin_is_frequency)
 See DissipatingZone::configure(). More...
 
virtual double modified_phase_lag (int orbital_frequency_multiplier, int spin_frequency_multiplier, double forcing_frequency, Dissipation::QuantityEntry entry, double &above_lock_value) const
 Should return the tidal phase lag times the love number for the given tidal term (or one of its derivatives). More...
 
virtual double love_coefficient (int, int, Dissipation::QuantityEntry) const
 Should return the corresponding component of the love coefficient (Lai 2012 Equation 24). More...
 
virtual CombinedStoppingConditionstopping_conditions (BinarySystem &system, bool primary, unsigned zone_index)
 Conditions detecting the next possible discontinuities in the evolution due to this zone. More...
 
virtual void change_e_order (unsigned new_e_order, BinarySystem &system, bool primary, unsigned zone_index)
 Changes the order of the eccentricity expansion performed. More...
 
virtual bool dissipative () const
 Return true iff disspation has been defined for the zone. More...
 
virtual bool can_lock () const
 Should return true iff the zone's dissipation is discontinuous at zero frequency. More...
 
 ~BrokenPowerlawPhaseLagZone ()
 Cleanup. More...
 
- Public Member Functions inherited from Evolve::DissipatingZone
void set_evolution_rates (double angular_momentum, double inclination, double periapsis)
 Set evolution rates for this zone's properties for storing in the eveloution. More...
 
double forcing_frequency (int orbital_frequency_multiplier, int spin_frequency_multiplier, double orbital_frequency) const
 The tidal forcing frequency for the given term and orbital frequency. More...
 
void set_reference_zone_angmom (double reference_angmom)
 Defines the angular momentum of the reference zone for single body evolution. More...
 
double periapsis_evolution (const Eigen::Vector3d &orbit_torque, const Eigen::Vector3d &zone_torque, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV, const Eigen::Vector3d &orbit_torque_deriv=Eigen::Vector3d(), const Eigen::Vector3d &zone_torque_deriv=Eigen::Vector3d())
 The rate at which the periapsis of the orbit/reference zone in this zone's equatorial plane is changing. More...
 
double inclination_evolution (const Eigen::Vector3d &orbit_torque, const Eigen::Vector3d &zone_torque, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV, const Eigen::Vector3d &orbit_torque_deriv=Eigen::Vector3d(), const Eigen::Vector3d &zone_torque_deriv=Eigen::Vector3d())
 The rate at which the inclination between this zone and the orbit is changing. More...
 
virtual bool locked (int orbital_frequency_multiplier, int spin_frequency_multiplier) const
 Should return true iff the given term is presently locked. More...
 
virtual bool locked () const
 Should return true iff any tidal term is locked. More...
 
const SpinOrbitLockInfolock_held () const
 The currntly held lock. More...
 
void release_lock ()
 Update the zone as necessary when the held lock disappears from the expansion. More...
 
void release_lock (short direction)
 Update the zone as necessary when the held lock is broken. More...
 
void set_lock (int orbital_frequency_multiplier, int spin_frequency_multiplier)
 Locks the zone spin to the orbit in the given ratio. More...
 
virtual double moment_of_inertia (int deriv_order=0) const =0
 Moment of inertia of the zone or its age derivative at the age of last configure() call. More...
 
virtual double moment_of_inertia (double age, int deriv_order=0) const =0
 The moment of inertia of the zone or its age derivative at a specified age (no configure necessary). More...
 
double spin_frequency () const
 The spin frequency of the given zone. More...
 
double angular_momentum () const
 The angular momentum of the given zone in \(M_\odot R_\odot^2\). More...
 
double tidal_power (bool above, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV) const
 The dimensionless tidal power or one of its derivatives. More...
 
double tidal_power (double above_fraction, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV) const
 Same as tidal_power(bool, Dissipation::QuantityEntry), but using the predefined mix of below/above contributions. More...
 
double tidal_torque_x (bool above, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV) const
 The dimensionless tidal torque along x. More...
 
double tidal_torque_x (double above_fraction, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV) const
 Same as tidal_torque_x(bool, Dissipation::QuantityEntry) but below and above contributions mixed. More...
 
double tidal_torque_y (bool above, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV) const
 The dimensionless torque along y. More...
 
double tidal_torque_y (double above_fraction, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV) const
 Same as tidal_torque_y(bool, Dissipation::QuantityEntry) but below and above contributions mixed. More...
 
double tidal_torque_z (bool above, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV) const
 The dimensionless tidal torque along z. More...
 
double tidal_torque_z (double above_fraction, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV) const
 Same as tidal_torque_z(bool, Dissipation::QuantityEntry) but below and above contributions mixed. More...
 
virtual double outer_radius (int deriv_order=0) const =0
 Outer radius of the zone or its derivative (per last. More...
 
virtual double outer_radius (double age, int deriv_order=0) const =0
 Same as outer_radius(int) but may be evaluated at a different age than for last confgure(). More...
 
virtual double outer_mass (int deriv_order=0) const =0
 Mass coordinate of the zone's outer ouboundary or its derivative. More...
 
virtual double outer_mass (double age, int deriv_order=0) const =0
 Same as outer_mass(int), but may be evaluated at a different age than last configure(). More...
 
virtual unsigned eccentricity_order () const
 
virtual void add_to_evolution ()
 Appends the state defined by last configure(), to the evolution. More...
 
virtual void rewind_evolution (unsigned nsteps)
 Discards the last steps from the evolution. More...
 
virtual void reset_evolution ()
 Discards all evolution. More...
 
const std::list< double > & get_evolution_real (ZoneEvolutionQuantities quantity) const
 The tabulated evolution of a real valued quantity so far. More...
 
const std::list< int > & get_evolution_integer (ZoneEvolutionQuantities quantity) const
 The tabulated evolution of an integer quantity so far. More...
 
unsigned locked_zone_index () const
 The index of this zone in the list of locked zones (valid only if locked). More...
 
unsigned & locked_zone_index ()
 Reference to the locked_zone_index() of this zone. More...
 
virtual void spin_jumped ()
 Notifies the zone that its spin just jumped discontinously. More...
 
virtual void reached_critical_age (double)
 Change the body as necessary at the given age. More...
 
virtual double next_stop_age () const
 The next age when the evolution needs to be stopped for a change in one of the bodies. More...
 
const SpinOrbitLockInfolock_monitored (bool other=false) const
 Useful for debugging. More...
 
- Public Member Functions inherited from Evolve::ZoneOrientation
 ZoneOrientation (double inclination=Core::NaN, double periapsis=Core::NaN)
 
void configure (double inclination, double periapsis)
 Changes the zone orientation. More...
 
void set_evolution_rates (double inclination, double periapsis)
 
double inclination (bool evolution_rate=false) const
 The angle between the angular momenta of the zone and the orbit. More...
 
double periapsis (bool evolution_rate=false) const
 The argument of periapsis of this zone minus the reference zone's. More...
 

Private Member Functions

std::vector< std::vector< double >::size_type >::size_type tidal_term_index (int orbital_frequency_multiplier, int spin_frequency_multiplier) const
 The index of the given tidal term within __tidal_indices. More...
 
double get_orbital_frequency (const BinarySystem &system) const
 Return the current orbital frequency of the given system. More...
 
double get_inertial_mode_factor (double abs_forcing_frequency, int orbital_frequency_multiplier, int spin_frequency_multiplier, Dissipation::QuantityEntry entry=Dissipation::NO_DERIV) const
 Calculate the factor by which dissipation is enhanced due to inertial modes or one of its logarithmic derivatives. More...
 
void reset ()
 
void set_spin_index ()
 Properly set the value of __spin_index per the current spin of the zone. More...
 
std::vector< double >::size_type get_tidal_index (double abs_forcing_frequency) const
 The index within __tidal_frequency_powers to use for the given forcing frequency. More...
 
void add_tidal_frequency_conditions (BinarySystem &system, bool primary, unsigned zone_index, CombinedStoppingCondition &result)
 Make sure that the entries in __tidal_frequency_conditions are appropriate for the current eccentricity expansion order. More...
 
void print_configuration (std::ostream &out_stream=std::clog)
 Print the configuration of the zone to stdlog. More...
 

Private Attributes

bool __dissipative
 Is the zone dissipative. More...
 
bool __can_lock
 
std::vector< double > __tidal_frequency_breaks
 The locations of the breaks in tidal frequency in rad/day. More...
 
std::vector< double > __spin_frequency_breaks
 The locations of the breaks in spin frequency in rad/day. More...
 
std::vector< double > __tidal_frequency_powers
 The powerlaw indices for the tidal frequency dependence. More...
 
std::vector< double > __spin_frequency_powers
 The powerlaw indices for the spin frequency dependence. More...
 
std::vector< double > __break_phase_lags
 The phase lags at the tidal/spin frequency breaks. More...
 
double __inertial_mode_enhancement
 See setup() More...
 
double __inertial_mode_sharpness
 See setup() More...
 
std::vector< double >::size_type __spin_index
 The index within __spin_frequency_powers of the powerlaw now in effect. More...
 
std::vector< std::vector< double >::size_type > __tidal_indices
 The indices within __tidal_frequency_powers of the powerlaw now in effect for each active tidal term. More...
 

Friends

class LagForcingFrequencyBreakCondition
 
class LagSpinBreakCondition
 

Additional Inherited Members

- Protected Member Functions inherited from Evolve::DissipatingZone
void initializing (bool flag)
 Notify the zone that it is in the process of initializing or not. More...
 
void configure_spin (double spin, bool spin_is_frequency)
 Configures only the spin of the zone. More...
 

Detailed Description

A DissipatingZone where the phase lag is described by a broken powerlaw.

By default the zone is non-dissipative.

Definition at line 26 of file BrokenPowerlawPhaseLagZone.h.

Constructor & Destructor Documentation

◆ BrokenPowerlawPhaseLagZone()

Evolve::BrokenPowerlawPhaseLagZone::BrokenPowerlawPhaseLagZone ( )
inline

Create a non-dissipative zone. Must call setup() if the zone is dissipative.

Definition at line 174 of file BrokenPowerlawPhaseLagZone.h.

◆ ~BrokenPowerlawPhaseLagZone()

Evolve::BrokenPowerlawPhaseLagZone::~BrokenPowerlawPhaseLagZone ( )
inline

Cleanup.

Definition at line 318 of file BrokenPowerlawPhaseLagZone.h.

Member Function Documentation

◆ add_tidal_frequency_conditions()

void Evolve::BrokenPowerlawPhaseLagZone::add_tidal_frequency_conditions ( BinarySystem system,
bool  primary,
unsigned  zone_index,
CombinedStoppingCondition result 
)
private

Make sure that the entries in __tidal_frequency_conditions are appropriate for the current eccentricity expansion order.

Parameters
systemThe system being evolved.
primaryIs the body this zone is part of, the primary in the system.
zone_indexThe index of the zone in the body.
resultThe final stopping condition returned by stopping_conditions() The newly constructed conditions get added to this argument.

Definition at line 131 of file BrokenPowerlawPhaseLagZone.cpp.

◆ can_lock()

virtual bool Evolve::BrokenPowerlawPhaseLagZone::can_lock ( ) const
inlinevirtual

Should return true iff the zone's dissipation is discontinuous at zero frequency.

Implements Evolve::DissipatingZone.

Definition at line 314 of file BrokenPowerlawPhaseLagZone.h.

◆ change_e_order()

void Evolve::BrokenPowerlawPhaseLagZone::change_e_order ( unsigned  new_e_order,
BinarySystem system,
bool  primary,
unsigned  zone_index 
)
virtual

Changes the order of the eccentricity expansion performed.

Parameters
new_e_orderThe new eccentricity expansion order.
systemThe system being evolved.
primaryIs the body this zone is part of, the primary in the system.
zone_indexThe index of the zone in the body.

Reimplemented from Evolve::DissipatingZone.

Definition at line 616 of file BrokenPowerlawPhaseLagZone.cpp.

◆ configure()

void Evolve::BrokenPowerlawPhaseLagZone::configure ( bool  initialize,
double  age,
double  orbital_frequency,
double  eccentricity,
double  orbital_angmom,
double  spin,
double  inclination,
double  periapsis,
bool  spin_is_frequency 
)
virtual

◆ dissipative()

virtual bool Evolve::BrokenPowerlawPhaseLagZone::dissipative ( ) const
inlinevirtual

Return true iff disspation has been defined for the zone.

Implements Evolve::DissipatingZone.

Definition at line 311 of file BrokenPowerlawPhaseLagZone.h.

◆ get_inertial_mode_factor()

double Evolve::BrokenPowerlawPhaseLagZone::get_inertial_mode_factor ( double  abs_forcing_frequency,
int  orbital_frequency_multiplier,
int  spin_frequency_multiplier,
Dissipation::QuantityEntry  entry = Dissipation::NO_DERIV 
) const
private

Calculate the factor by which dissipation is enhanced due to inertial modes or one of its logarithmic derivatives.

Parameters
abs_forcing_frequencyAbsolute value of the frequency of the tidal term for which
orbital_frequency_multiplierThe multiplier of the orbital frequency in the expression for the forcing frequency. Only used if derivative is requested.
spin_frequency_multiplierThe multiplier of the spin frequency in the expression for the forcing frequency. Only used if derivative is requested.
entryIf Dissipation::NO_DERIV, return the enhancement factor. Ifthe specified quantity. If anything else, derivative is zero.

Definition at line 18 of file BrokenPowerlawPhaseLagZone.cpp.

◆ get_orbital_frequency()

double Evolve::BrokenPowerlawPhaseLagZone::get_orbital_frequency ( const BinarySystem system) const
private

Return the current orbital frequency of the given system.

Definition at line 7 of file BrokenPowerlawPhaseLagZone.cpp.

◆ get_tidal_index()

std::vector< double >::size_type Evolve::BrokenPowerlawPhaseLagZone::get_tidal_index ( double  abs_forcing_frequency) const
private

The index within __tidal_frequency_powers to use for the given forcing frequency.

Should only be called when starting the evolution in a two-body configuration. After that, this values should be updated by stopping conditions.

Parameters
abs_forcing_frequencyThe absolute value of the forcing frequency.

Definition at line 112 of file BrokenPowerlawPhaseLagZone.cpp.

◆ love_coefficient()

virtual double Evolve::BrokenPowerlawPhaseLagZone::love_coefficient ( int  ,
int  ,
Dissipation::QuantityEntry   
) const
inlinevirtual

Should return the corresponding component of the love coefficient (Lai 2012 Equation 24).

Implements Evolve::DissipatingZone.

Definition at line 265 of file BrokenPowerlawPhaseLagZone.h.

◆ modified_phase_lag()

double Evolve::BrokenPowerlawPhaseLagZone::modified_phase_lag ( int  orbital_frequency_multiplier,
int  spin_frequency_multiplier,
double  forcing_frequency,
Dissipation::QuantityEntry  entry,
double &  above_lock_value 
) const
virtual

Should return the tidal phase lag times the love number for the given tidal term (or one of its derivatives).

In case the forcing frequency is exactly zero, it should return the phase lag for the case of the spin frequency approaching the term from below. The lag for spin frequency approaching from above should be written to above_lock_value. If the forcing frequency is non-zero, leave above_lock_value untouched.

Parameters
orbital_frequency_multiplierThe multiplier of the orbital frequency in the expression for the forcing frequency.
spin_frequency_multiplierThe multiplier of the spin frequency in the expression for the forcing frequency.
forcing_frequencyThe current forcing frequency in rad/day.
entryThe return value should be either the phase lag itself (NO_DERIV) or its derivative w.r.t. the specified quantity.
above_lock_valueIf the lag of a locked term is calculated this should be set to the lag assuming the spin frequency is just above the lock. Otherwise, leave untouched.

Implements Evolve::DissipatingZone.

Definition at line 419 of file BrokenPowerlawPhaseLagZone.cpp.

◆ print_configuration()

void Evolve::BrokenPowerlawPhaseLagZone::print_configuration ( std::ostream &  out_stream = std::clog)
private

Print the configuration of the zone to stdlog.

Definition at line 164 of file BrokenPowerlawPhaseLagZone.cpp.

◆ reset()

void Evolve::BrokenPowerlawPhaseLagZone::reset ( )
private

Cleanup either during destruction or in preparation for changing dissipation.

Definition at line 74 of file BrokenPowerlawPhaseLagZone.cpp.

◆ set_spin_index()

void Evolve::BrokenPowerlawPhaseLagZone::set_spin_index ( )
private

Properly set the value of __spin_index per the current spin of the zone.

Should only be called once, when starting the evolution in a two-body configuration. After that, this values should be updated by stopping conditions.

Definition at line 83 of file BrokenPowerlawPhaseLagZone.cpp.

◆ setup()

void Evolve::BrokenPowerlawPhaseLagZone::setup ( const std::vector< double > &  tidal_frequency_breaks,
const std::vector< double > &  spin_frequency_breaks,
const std::vector< double > &  tidal_frequency_powers,
const std::vector< double > &  spin_frequency_powers,
double  reference_phase_lag,
double  inertial_mode_enhancement = 1.0,
double  inertial_mode_sharpness = 10.0 
)

Seup the zone with the given breaks/powers and inertial mode enhancement. Continuous accress all breaks. Must only be called before use.

Parameters
tidal_frequency_breaksThe locations of the breaks in tidal frequency in rad/day. Entries should be sorted.
spin_frequency_breaksThe locations of the breaks in spin frequency in rad/day. Entries should be sorted.
tidal_frequency_powersThe powerlaw indices for the tidal frequency dependence. Should be indexed in the same order as tidal_frequency_breaks, but must contain an additional starting entry for the powerlaw index before the first break.
spin_frequency_powersThe powerlaw indices for the spin frequency dependence. Should be indexed in the same order as spin_frequency_breaks, but must contain an additional starting entry for the powerlaw index before the first break.
reference_phase_lagThe phase lag at the first tidal and first spin frequency break. The rest are calculated by imposing continuity.
inertial_mode_enhancementThe factor by which dissipation is enhanced in the inertial mode range.
inertial_mode_sharpnessA parameter controlling how sharp the transition between non-enhanced and enhanced dissipation is near the inertial mode boundary.

The expression for the enhancement factor is: \( \max\left[\left|\frac{2\Omega_\star}{\Omega_{m,s}}\right|^\beta, \gamma\right] \) Where \(\beta\) is inertial_mode_sharpness, and \(\gamma\) is inertial_mode_enhancement.

Definition at line 218 of file BrokenPowerlawPhaseLagZone.cpp.

◆ stopping_conditions()

CombinedStoppingCondition * Evolve::BrokenPowerlawPhaseLagZone::stopping_conditions ( BinarySystem system,
bool  primary,
unsigned  zone_index 
)
virtual

Conditions detecting the next possible discontinuities in the evolution due to this zone.

Must be deleted when no longer necessary.

Parameters
systemThe system being evolved.
primaryIs the body this zone is part of, the primary in the system.
zone_indexThe index of the zone in the body.

Reimplemented from Evolve::DissipatingZone.

Definition at line 583 of file BrokenPowerlawPhaseLagZone.cpp.

◆ tidal_term_index()

std::vector< std::vector<double>::size_type >::size_type Evolve::BrokenPowerlawPhaseLagZone::tidal_term_index ( int  orbital_frequency_multiplier,
int  spin_frequency_multiplier 
) const
inlineprivate

The index of the given tidal term within __tidal_indices.

Parameters
orbital_frequency_multiplierThe orbital frequency multiplier of the forcing term being monitored.
spin_frequency_multiplierThe spin frequency multiplier of the forcing term being monitored.

Definition at line 80 of file BrokenPowerlawPhaseLagZone.h.

Member Data Documentation

◆ __break_phase_lags

std::vector<double> Evolve::BrokenPowerlawPhaseLagZone::__break_phase_lags
private

The phase lags at the tidal/spin frequency breaks.

The tidal frequency break index changes faster and the spin frequency breaks index changes slower.

Definition at line 40 of file BrokenPowerlawPhaseLagZone.h.

◆ __can_lock

bool Evolve::BrokenPowerlawPhaseLagZone::__can_lock
private

Is the zone's dissipation discontinuous at zero frequency (hence spin-orbit locks are a possibility).

Definition at line 32 of file BrokenPowerlawPhaseLagZone.h.

◆ __dissipative

bool Evolve::BrokenPowerlawPhaseLagZone::__dissipative
private

Is the zone dissipative.

Definition at line 32 of file BrokenPowerlawPhaseLagZone.h.

◆ __inertial_mode_enhancement

double Evolve::BrokenPowerlawPhaseLagZone::__inertial_mode_enhancement
private

See setup()

Definition at line 59 of file BrokenPowerlawPhaseLagZone.h.

◆ __inertial_mode_sharpness

double Evolve::BrokenPowerlawPhaseLagZone::__inertial_mode_sharpness
private

See setup()

Definition at line 59 of file BrokenPowerlawPhaseLagZone.h.

◆ __spin_frequency_breaks

std::vector<double> Evolve::BrokenPowerlawPhaseLagZone::__spin_frequency_breaks
private

The locations of the breaks in spin frequency in rad/day.

Definition at line 40 of file BrokenPowerlawPhaseLagZone.h.

◆ __spin_frequency_powers

std::vector<double> Evolve::BrokenPowerlawPhaseLagZone::__spin_frequency_powers
private

The powerlaw indices for the spin frequency dependence.

Definition at line 40 of file BrokenPowerlawPhaseLagZone.h.

◆ __spin_index

std::vector<double>::size_type Evolve::BrokenPowerlawPhaseLagZone::__spin_index
private

The index within __spin_frequency_powers of the powerlaw now in effect.

Definition at line 66 of file BrokenPowerlawPhaseLagZone.h.

◆ __tidal_frequency_breaks

std::vector<double> Evolve::BrokenPowerlawPhaseLagZone::__tidal_frequency_breaks
private

The locations of the breaks in tidal frequency in rad/day.

Definition at line 40 of file BrokenPowerlawPhaseLagZone.h.

◆ __tidal_frequency_powers

std::vector<double> Evolve::BrokenPowerlawPhaseLagZone::__tidal_frequency_powers
private

The powerlaw indices for the tidal frequency dependence.

Definition at line 40 of file BrokenPowerlawPhaseLagZone.h.

◆ __tidal_indices

std::vector< std::vector<double>::size_type > Evolve::BrokenPowerlawPhaseLagZone::__tidal_indices
private

The indices within __tidal_frequency_powers of the powerlaw now in effect for each active tidal term.

The index of the (m' * Worb - m * Wspin) term is given by:

  • (2 + m) + 5 * m' if m' > 0
  • (2 - m) - 5 * m' if m' < 0 This way changing the eccentricity expansion coefficient only affects the tail of the vector.

Definition at line 76 of file BrokenPowerlawPhaseLagZone.h.


The documentation for this class was generated from the following files: