  We report on the calculation of the threshold soft function for heavy quark
pair production in e+ e- annihilation at two-loop order. Our main result is a
generalization of the familiar Drell-Yan threshold soft function to the case of
non-zero primary quark mass. We set up a framework based on the method of
differential equations which allows for the straightforward calculation of the
bare soft function to arbitrarily high orders in the dimensional regularization
parameter. Remarkably, we find that we can obtain the bare two-loop Drell-Yan
soft function from the heavy quark soft function to the order in epsilon
required for a two-loop calculation by making simple replacements. We expect
that our results will be of use, both as an important input for precision
physics calculations at linear colliders and, more formally, as a first step
towards a better understanding of the connection between vacuum matrix elements
of massive soft Wilson lines and vacuum matrix elements of massless soft Wilson
lines.
