public static function Ld($bm, array $p, array $q, array $e, $em, $dlim,
          array &$p1) {

    $i;
    $qpe = [];
    $qdqpe;
    $w;
    $eq  = [];
    $peq = [];

    /* q . (q + e). */
    for ($i = 0; $i < 3; $i++) {
      $qpe[$i] = $q[$i] + $e[$i];
    }
    $qdqpe = IAU::Pdp($q, $qpe);

    /* 2 x G x bm / ( em x c^2 x ( q . (q + e) ) ). */
    $w = $bm * SRS / $em / max($qdqpe, $dlim);

    /* p x (e x q). */
    IAU::Pxp($e, $q, $eq);
    IAU::Pxp($p, $eq, $peq);

    /* Apply the deflection. */
    for ($i = 0; $i < 3; $i++) {
      $p1[$i] = $p[$i] + $w * $peq[$i];
    }

    /* Finished. */
  }