docs: def withinroche(semimajors,M1,R1,M2,R2):
      """
      Returns boolean array that is True where two stars are within Roche lobe
      """
      q = M1/M2
      return ((R1+R2)*RSUN) > (rochelobe(q)*semimajors*AU)