Constraints

Setup

Let's first make sure we have the latest version of PHOEBE 2.3 installed (uncomment this line if running in an online notebook session such as colab).


In [1]:
#!pip install -I "phoebe>=2.3,<2.4"

In [2]:
import phoebe
from phoebe import u # units
import numpy as np
import matplotlib.pyplot as plt

logger = phoebe.logger()

b = phoebe.default_binary()

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 by filtering on context='constraint'.


In [3]:
b.filter(context='constraint')


Out[3]:
<ParameterSet: 26 parameters | kinds: extinction, star, orbit>

To see what all of these constraints do, see Advanced: Built-In Constraints or look at the constraint API docs.

For now let's look at a single constraint by accessing a ConstraintParameter.


In [4]:
b['constraint']['primary']['mass']


Out[4]:
<ConstraintParameter: {mass@primary@component} = (39.478418 * ({sma@binary@component} ** 3.000000)) / ((({period@binary@component} ** 2.000000) * (1.000000 + {q@binary@component})) * 2942.206217504418873431859537959099) (solar units) => 0.9988131358058302 solMass>

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.


In [5]:
print(b.get_value('mass@primary@component'))


0.9988131358058302

The parameter is aware that it's being constrained and has references to all the relevant linking parameters.


In [6]:
print(b['mass@primary@component'])


Parameter: mass@primary@component
                       Qualifier: mass
                     Description: Mass
                           Value: 0.9988131358058302 solMass
                  Constrained by: sma@binary@component, period@binary@component, q@binary@component
                      Constrains: logg@primary@component, mass@secondary@component
                      Related to: requiv@primary@component, logg@primary@component, sma@binary@component, period@binary@component, q@binary@component, mass@secondary@component

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 Advanced: Constraints and Changing Hierarchices for more details.

Re-Parameterizing or "Flipping" Constraints

NOTE: 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 we would rather provide masses for some reason (perhaps that was what was provided in a paper). We can choose to provide mass and instead have one of its related parameters constrained by calling b.flip_constraint.


In [7]:
print(b['mass@primary@component'].constrained_by)


[<Parameter: sma=5.3 solRad | keys: description, value, quantity, default_unit, limits, visible_if, copy_for, readonly, advanced>, <Parameter: period=1.0 d | keys: description, value, quantity, default_unit, limits, visible_if, copy_for, readonly, advanced>, <Parameter: q=1.0 | keys: description, value, quantity, default_unit, limits, visible_if, copy_for, readonly, advanced>]

In [8]:
print(b['value@mass@primary@component'], b['value@mass@secondary@component'], b['value@period@orbit@component'])


0.9988131358058302 0.9988131358058302 1.0

In [9]:
b.flip_constraint('mass@primary', solve_for='period')


Out[9]:
<ConstraintParameter: {period@binary@component} = ((39.478418 * ({sma@binary@component} ** 3.000000)) / (({mass@primary@component} * (1.000000 + {q@binary@component})) * 2942.206217504418873431859537959099)) ** (1./2) (solar units) => 1.0 d>

In [10]:
b['mass@primary@component'] = 1.0

In [11]:
print(b['value@mass@primary@component'], b['value@mass@secondary@component'], b['value@period@orbit@component'])


1.0 1.0000000000000002 0.9994063917175185

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.


In [12]:
print(b['constraint'])


ParameterSet: 26 parameters
                   ebv@constraint: {Av@system} / {Rv@system}
          freq@primary@constraint: 6.283185 / {period@primary@component}
          logg@primary@constraint: log10((({mass@primary@component} / ({requiv@primary@component} ** 2.000000)) * 0.000000000066740800000000003352) * 100.000000)
   irrad_frac_lost_bol@primary...: 1.000000 - {irrad_frac_refl_bol@primary@component}
        freq@secondary@constraint: 6.283185 / {period@secondary@component}
        logg@secondary@constraint: log10((({mass@secondary@component} / ({requiv@secondary@component} ** 2.000000)) * 0.000000000066740800000000003352) * 100.000000)
   irrad_frac_lost_bol@seconda...: 1.000000 - {irrad_frac_refl_bol@secondary@component}
                 asini@constraint: {sma@binary@component} * (sin({incl@binary@component}))
            t0_perpass@constraint: t0_supconj_to_perpass({t0_supconj@binary@component}, {period@binary@component}, {ecc@binary@component}, {per0@binary@component}, {dpdt@binary@component}, {dperdt@binary@component}, {t0@system})
                t0_ref@constraint: t0_supconj_to_ref({t0_supconj@binary@component}, {period@binary@component}, {ecc@binary@component}, {per0@binary@component}, {dpdt@binary@component}, {dperdt@binary@component}, {t0@system})
             mean_anom@constraint: (6.283185 * ({t0@system} - {t0_perpass@binary@component})) / {period@binary@component}
                 ecosw@constraint: {ecc@binary@component} * (cos({per0@binary@component}))
                 esinw@constraint: {ecc@binary@component} * (sin({per0@binary@component}))
           freq@binary@constraint: 6.283185 / {period@binary@component}
         period@binary@constraint: ((39.478418 * ({sma@binary@component} ** 3.000000)) / (({mass@primary@component} * (1.000000 + {q@binary@component})) * 2942.206217504418873431859537959099)) ** (1./2)
           sma@primary@constraint: {sma@binary@component} / (1.000000 + (1.000000 / {q@binary@component}))
        period@primary@constraint: {period@binary@component} / {syncpar@primary@component}
          incl@primary@constraint: {incl@binary@component} + {pitch@primary@component}
       long_an@primary@constraint: {long_an@binary@component} + {yaw@primary@component}
    requiv_max@primary@constraint: requiv_L1({q@binary@component}, {syncpar@primary@component}, {ecc@binary@component}, {sma@binary@component}, {incl@primary@component}, {long_an@primary@component}, {incl@binary@component}, {long_an@binary@component}, 1)
                  mass@constraint: (39.478418 * ({sma@binary@component} ** 3.000000)) / ((({period@binary@component} ** 2.000000) * (1.000000 + (1.000000 / {q@binary@component}))) * 2942.206217504418873431859537959099)
         sma@secondary@constraint: {sma@binary@component} / (1.000000 + {q@binary@component})
      period@secondary@constraint: {period@binary@component} / {syncpar@secondary@component}
        incl@secondary@constraint: {incl@binary@component} + {pitch@secondary@component}
     long_an@secondary@constraint: {long_an@binary@component} + {yaw@secondary@component}
   requiv_max@secondary@constr...: requiv_L1({q@binary@component}, {syncpar@secondary@component}, {ecc@binary@component}, {sma@binary@component}, {incl@secondary@component}, {long_an@secondary@component}, {incl@binary@component}, {long_an@binary@component}, 2)

In [13]:
b['period@constraint@binary']


Out[13]:
<ConstraintParameter: {period@binary@component} = ((39.478418 * ({sma@binary@component} ** 3.000000)) / (({mass@primary@component} * (1.000000 + {q@binary@component})) * 2942.206217504418873431859537959099)) ** (1./2) (solar units) => 0.9994063917175185 d>

In [14]:
b['period@constraint@binary'].meta


Out[14]:
OrderedDict([('time', None),
             ('qualifier', 'period'),
             ('history', None),
             ('feature', None),
             ('component', 'binary'),
             ('dataset', None),
             ('constraint', None),
             ('distribution', None),
             ('compute', None),
             ('model', None),
             ('solver', None),
             ('solution', None),
             ('figure', None),
             ('kind', 'orbit'),
             ('context', 'constraint'),
             ('twig', 'period@binary@orbit@constraint'),
             ('uniquetwig', 'period@binary@constraint')])

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).

Next

Next up: let's add a dataset to our Bundle.

Or look at any of the advanced constraints topics: