def keplers_third_law_hierarchical(b, orbit1, orbit2, solve_for=None, **kwargs):
    """
    TODO: add documentation

    This is usually passed as an argument to
     <phoebe.frontend.bundle.Bundle.add_constraint>.
    """

    hier = b.hierarchy

    orbit1_ps = _get_system_ps(b, orbit1)
    orbit2_ps = _get_system_ps(b, orbit2)

    sma1 = orbit1_ps.get_parameter(qualifier='sma', **_skip_filter_checks)
    sma2 = orbit2_ps.get_parameter(qualifier='sma', **_skip_filter_checks)

    q1 = orbit1_ps.get_parameter(qualifier='q', **_skip_filter_checks)
    q2 = orbit2_ps.get_parameter(qualifier='q', **_skip_filter_checks)

    period1 = orbit1_ps.get_parameter(qualifier='period', **_skip_filter_checks)
    period2 = orbit2_ps.get_parameter(qualifier='period', **_skip_filter_checks)

    # NOTE: orbit1 is the outer, so we need to check orbit2... which will
    # be the OPPOSITE component as that of the mass we're solving for
    if hier.get_primary_or_secondary(orbit2_ps.component) == 'primary':
        # NOTE: 1.0 must be on the right for distribution math to behave
        qthing1 = q1+1.0
    else:
        # NOTE: 1.0 must be on the right for distribution math to behave
        qthing1 = 1./q1+1.0

    if solve_for in [None, sma1]:
        lhs = sma1
        rhs = (sma2**3 * qthing1 * period1**2/period2**2)**"(1./3)"
    else:
        # TODO: add other options to solve_for
        raise NotImplementedError

    return lhs, rhs, [], {'orbit1': orbit1, 'orbit2': orbit2}