def keplersconstant(m1, m2):
    """Compute Kepler's constant for two gravitationally bound masses
        Inputs:
            m1,m2 (arr,flt): masses of the two objects in solar masses.  Must be astropy objects
        Returns:
            Kepler's constant in m^3 s^(-2)
    """
    import astropy.constants as c
    import astropy.units as u
    m1 = m1.to(u.Msun)
    m2 = m2.to(u.Msun)
    mu = c.G * m1 * m2
    m = (1 / m1 + 1 / m2) ** -1
    kep = mu / m
    return kep.to(u.m ** 3 / u.s ** 2)
