At the very least, you could do it with logs like so:

    >>> from math import exp, log
    >>> ref  = [i for i in range(1, 500000)]
    >>> test = [i if i != 7922 else 59148 for i in ref]
    >>> s = sum(test) - sum(ref)
    >>> r = exp(sum(map(log, test)) - sum(map(log, ref)))
    >>> b = s/(r - 1); b
    7921.999994157095
    >>> (int(round(b)), s + int(round(b)))
    (7922, 59148)

The use of logs is really just a convenience for finding product(test) and product(ref), so that we have both a - b (from the sums) and a/b (from the ratio of the products).

Python of course has bignums, so we can also calculate these products directly, but that gets to be a very large calculation; you don't want to be storing N! for N ~= 2^64 if at all possible.

If we know that these numbers are bounded by 2^n (e.g. they're all longs) then we might be able to do this product modulo some larger base, perhaps even 2^n, with the assurance that when we do the division all of that crap will cancel; I don't know and haven't tried it.