!pip install -I "phoebe>=2.0,<2.1"

"""
Explanation: Constraints
Setup
Let's first make sure we have the latest version of PHOEBE 2.0 installed. (You can comment out this line if you don't use pip for your installation or don't want to update to the latest release).
End of explanation
"""


import phoebe
from phoebe import u # units
import numpy as np
import matplotlib.pyplot as plt

logger = phoebe.logger()

b = phoebe.default_binary()

"""
Explanation: As always, let's do imports and initialize a logger and a new Bundle.  See Building a System for more details.
End of explanation
"""


b.filter(context='constraint')

"""
Explanation: What are Constraints?
Constraints live in their own context of the Bundle, and many are created
by default - either when you add a component or when you set the system hierarchy.
Let's look at all the existing constraints for our binary system.
End of explanation
"""


b['constraint']['primary']['mass']

"""
Explanation: To see what all of these constraints do, see the 'Built-in Constraints' section below.
For now let's look at a single constraint.
End of explanation
"""


print b.get_value('mass@primary@component')

"""
Explanation: Here we see the equation used to derive the mass of the primary star
from its orbit, as well as the current value
If we access the Parameter that it is constraining we can see that it
is automatically kept up-to-date.
End of explanation
"""


print b['mass@primary@component']

"""
Explanation: The parameter is aware that it's being
constrained and all the relevant linking parameters.
End of explanation
"""


b['asini@constraint']

"""
Explanation: If you change the hierarchy, built-in cross-object constraints (like mass
that depends on its parent orbit), will be adjusted to reflect the new hierarchy.  See the 'Changing Hierarchies' section below for more details.
Built-in Constraints
There are a number of built-in constraints that will be applied to your system by default.  These are all listed below:
asini
This constraint handles computing the projected semi-major axis along the line of sight and can be automatically inverted to solve for either 'asini', 'sma', or 'incl'.
End of explanation
"""


b['esinw@constraint']

b['ecosw@constraint']

"""
Explanation: esinw, ecosw
These constraints handle computing the projected eccentricity which can be helpful in that they are better representations of the geometry of a light curve and result in symmetric posteriors for near-circular orbits.
Both can be inverted to also automatically solve for 'ecc' or 'per0'.
End of explanation
"""


b['t0_perpass@constraint']

"""
Explanation: t0
This constraint handles converting between different t0 conventions - namely providing a reference time at periastron passage (t0_perpass) and at superior conjunction (t0_supconj).
Currently, this constraint only supports inverting to be solved for 't0_supconj' (ie you cannot automatically invert this constraint to constraint phshift or per0).
End of explanation
"""


b['freq@constraint']

b['freq@binary@constraint']

b['freq@primary@constraint']

"""
Explanation: freq
This constraint handles the simple conversion to frequency from period - whether that be rotational or orbital - and does support inversion to solve for 'period'.
End of explanation
"""


b['mass@constraint']

b['mass@primary@constraint']

"""
Explanation: mass
This constraint handles solving for the mass of a component by obeying Kepler's third law within the parent orbit.
It can be inverted to solve for 'sma' or 'period' (in addition to 'mass'), but not 'q'.
End of explanation
"""


b['sma@constraint']

b['sma@primary@constraint']

"""
Explanation: component sma
This constraint handles computing the semi-major axis of a component about the center of mass of its parent orbit.  Note that this is not the same as the semi-major axis of the parent orbit.
This currently can be inverted to solve for 'sma' of the parent orbit, but not 'q'.
End of explanation
"""


b['pot@constraint']

b['pot@primary@constraint']

"""
Explanation: potential
This constraint handles solving for the equipotential of the roche surface from the polar radius.
This currently can be inverted to solve for 'rpole' but not 'sma', 'syncpar', 'ecc', or 'q'.
End of explanation
"""


b['period@constraint']

b['period@primary@constraint']

"""
Explanation: rotation period
This constraint handles computing the rotation period of a star given its synchronicity parameter (syncpar).
It can be inverted to solve for any of the three parameters 'period' (both rotational and orbital) and 'syncpar'.
End of explanation
"""


b['incl@constraint']

b['incl@primary@constraint']

"""
Explanation: inclination (aligned binary)
This constraint handles requiring the inclination of a star to be the same as that of the parent orbit.  With the constraint enabled, the binary is required to be "aligned".  Disabling this constraint (which currently is not yet supported) will allow for misaligned binaries.
This constraint can be inverted to solve for either 'incl'.
End of explanation
"""


print b['mass@primary@component'].constrained_by

print b['value@mass@primary@component'], b['value@mass@secondary@component'], b['value@period@orbit@component']

b.flip_constraint('mass@primary', 'period')

b['mass@primary@component'] = 1.0

print b['value@mass@primary@component'], b['value@mass@secondary@component'], b['value@period@orbit@component']

"""
Explanation: Re-Parameterizing
NOTE: this is an experimental feature.  When re-parameterizing, please be
careful and make sure all results and parameters make sense.
As we've just seen, the mass is a constrained (ie derived) parameter.  But
let's say that you would rather provide masses for some reason (perhaps
that was what was provided in a paper).  You can choose to provide mass
and instead have one of its related parameters constrained
End of explanation
"""


print b['constraint']

b['period@constraint@binary']

b['period@constraint@binary'].meta

"""
Explanation: You'll see that when we set the primary mass, the secondary mass has also changed (because the masses are related through q) and the period has changed (based on resolving the Kepler's third law constraint).
Note that the tags for the constraint are based on those of the constrained parameter, so to switch the parameterization back, we'll have to use a slightly different twig.
End of explanation
"""


b.flip_constraint('period@binary', 'mass')

"""
Explanation: Notice that the qualifier tag has changed from 'mass' to 'period' and the 'component' tag has changed from 'primary' to 'binary' (since sma is in the binary).
End of explanation
"""


b.set_value('q', 0.8)

"""
Explanation: Changing Hierarchies
NOTE: this is an experimental feature.  Just as when re-parameterizing, please be
careful and make sure all results and parameters make sense.
Some of the built-in constraints depend on the system hierarchy, and will automatically adjust to reflect changes to the hierarchy.
For example, the masses depend on the period and semi-major axis of the parent orbit but also depend on the mass-ratio (q) which is defined as the primary mass over secondary mass.  For this reason, changing the roles of the primary and secondary components should be reflected in the masses (so long as q remains fixed).
In order to show this example, let's set the mass-ratio to be non-unity.
End of explanation
"""


print "M1: {}, M2: {}".format(b.get_value('mass@primary@component'),
                              b.get_value('mass@secondary@component'))

"""
Explanation: Here the star with component tag 'primary' is actually the primary component in the hierarchy, so should have the LARGER mass (for a q < 1.0).
End of explanation
"""


b.set_hierarchy('orbit:binary(star:secondary, star:primary)')

print b.get_value('q')

print "M1: {}, M2: {}".format(b.get_value('mass@primary@component'),
                              b.get_value('mass@secondary@component'))

"""
Explanation: Now let's flip the hierarchy so that the star with the 'primary' component tag is actually the secondary component in the system (and so takes the role of numerator in q = M2/M1).
For more information on the syntax for setting hierarchies, see the Building a System Tutorial.
End of explanation
"""


print "M1: {}, M2: {}, period: {}, q: {}".format(b.get_value('mass@primary@component'),
                                                 b.get_value('mass@secondary@component'),
                                                 b.get_value('period@binary@component'),
                                                 b.get_value('q@binary@component'))

b.flip_constraint('mass@secondary@constraint', 'period')

print "M1: {}, M2: {}, period: {}, q: {}".format(b.get_value('mass@primary@component'),
                                                 b.get_value('mass@secondary@component'),
                                                 b.get_value('period@binary@component'),
                                                 b.get_value('q@binary@component'))

b.set_value('mass@secondary@component', 1.0)

print "M1: {}, M2: {}, period: {}, q: {}".format(b.get_value('mass@primary@component'),
                                                 b.get_value('mass@secondary@component'),
                                                 b.get_value('period@binary@component'),
                                                 b.get_value('q@binary@component'))

"""
Explanation: Even though under-the-hood the constraints are being rebuilt from scratch, they will remember if you have flipped them to solve for some other parameter.
To show this, let's flip the constraint for the secondary mass to solve for 'period' and then change the hierarchy back to its original value.
End of explanation
"""


