Perhaps someone can verify the result below, or suggest cleaner ways to expand/write this?
Starting from the definition of a hyperbola as the locus of points X = (x,y) that have a constant distance difference to two points (the foci, A = (a_x, a_y) and B = (b_x, b_y)): $ |A-X| - |B-X| = d $
We can write this (at least half of the solution) as:
$ \sqrt{(a_x-x)^2 + (a_y-y)^2} - \sqrt{(b_x-x)^2 + (b_y-y)^2} = d $
There is tremendous manipulation necessary to write this equation in general quadratic form ($Ax^2 + By^2 + Cxy + Dx + Ey + F = 0$):
Move one of the square roots to the other side of the equation so that when we square both sides we will not eventually introduce quartic terms: $ \sqrt{(a_x-x)^2 + (a_y-y)^2} = \sqrt{(b_x-x)^2 + (b_y-y)^2} + d $
Square both sides: $ (a_x-x)^2 + (a_y-y)^2 = (b_x-x)^2 + (b_y-y)^2 + 2d \sqrt{(b_x-x)^2 + (b_y-y)^2} + d^2 $
Expand: $ a_x^2 -2xa_x + x^2 + a_y^2 - 2a_y y + y^2 = b_x^2 - 2xb_x + x^2 + b_y^2 -2b_y y + y^2 + 2d \sqrt{(b_x-x)^2 + (b_y-y)^2} + d^2 $
Cancel $x$ and $y$ quadratic terms (that occur on both sides of the equation): $ a_x^2 -2xa_x + a_y^2 - 2a_y y = b_x^2 - 2xb_x + b_y^2 -2b_y y + 2d \sqrt{(b_x-x)^2 + (b_y-y)^2} + d^2 $
Move everything not connect to the square root to the left side so we can eventually square both sides again to get rid of the square root: $ a_x^2 - 2xa_x + a_y^2 - 2a_y y - b_x^2 = 2xb_x - b_y^2 + 2b_y y - d^2 = 2d \sqrt{(b_x-x)^2 + (b_y-y)^2} $
Group terms on the left side into the form Jx + Ky + L to make squaring easier: $ \left((-2a_x + 2b_x)x + (-2a_y + 2b_y)y + (a_x^2 + a_y^2 - b_x^2 - b_y^2 - d^2)\right) = 2d \sqrt{(b_x-x)^2 + (b_y-y)^2} $
Assign $J = -2a_x + 2b_x$, $k = -2a_y + 2b_y$, $L = a_x^2 + a_y^2 - b_x^2 - b_y^2 - d^2$ so this can be written as
$ Jx + Ky + L = 2d \sqrt{(b_x-x)^2 + (b_y-y)^2} $
Square both sides:
$ (Jx + Ky + L)^2 = 4d^2 \left((b_x-x)^2 + (b_y-y)^2\right) $
$ J^2x^2 + JxKy + JxL + KyJx + K^2y^2 + KyL + LJx + LKy + L^2 = 4d^2 \left(b_x^2 - 2xb_x +x^2 + b_y^2 - 2yb_y + y^2\right) $
Combine terms on the left, and multiply through on the right:
$ J^2x^2 + K^2y^2 + 2JxKy + 2JxL + 2KyL + L^2 = 4d^2 b_x^2 - 8d^2xb_x + 4d^2x^2 + 4d^2b_y^2 - 8d^2yb_y + 4d^2y^2 $
Combine all terms ($x^2, y^2, xy, x, y$):
$ (J^2-4d^2)x^2 + (K^2-4d^2)y^2 + (2JK)xy + (2JL+8d^2b_x)x + (2KL + 8d^2b_y)y + (L^2 - 4d^2 b_x^2 - 4d^2 b_y^2) = 0 $
Substituting back in the values of J, K, and L:
$ ((-2a_x + 2b_x)^2-4d^2)x^2 + ((-2a_y + 2b_y)^2-4d^2)y^2 + (2(-2a_x + 2b_x)(-2a_y + 2b_y))xy \\+ (2(-2a_x + 2b_x)(a_x^2 + a_y^2 - b_x^2 - b_y^2 - d^2)+8d^2b_x)x + (2(-2a_y + 2b_y)(a_x^2 + a_y^2 - b_x^2 - b_y^2 - d^2) + 8d^2b_y)y \\+ ((a_x^2 + a_y^2 - b_x^2 - b_y^2 - d^2)^2 - 4d^2 b_x^2 - 4d^2 b_y^2) = 0 $
Expand: $ (4a_x^2 - 8a_x b_x + 4b_x^2)x^2 + \\ (4a_y^2 - 8a_y b_y + 4b_y^2)y^2 + \\ (8a_xa_y - 8a_xb_y - 8a_yb_x + 8b_xb_y)xy +\\ (-4a_x^3 -4a_xa_y^2 + 4 a_x b_x^2 + 4 a_x b_y^2 + 4 a_x d^2 + 4b_xa_x^2 + 4b_xa_y^2 + 4b_x^3 - 4b_xb_y^2 + 4b_xd^2)x +\\ (-4a_ya_x^2 - 4a_y^3 + 4 a_y b_x^2 + 4 a_y b_y^2 + 4 a_y d^2 + 4b_ya_x^2 + 4b_ya_y^2 - 4b_yb_x^2 - 4b_y^3 + 4b_yd^2)y +\\ a_x^4 + a_x^2 a_y^2 - a_x^2b_x^2 - a_x^2b_y^2 - a_x^2d^2\\ + a_y^2 a_x^2 + a_y^4 - a_y^2b_x^2 - a_y^2b_y^2 - a_y^2 d^2\\ - b_x^2 a_x^2 - b_x^2 a_y^2 + b_x^4 + b_x^2 b_y^2 -3 b_x^2 d_2\\ - b_y^2 a_x^2 - b_y^2a_y^2 + b_y^2b_x^2 + b_y^4 + b_y^2 d^2\\ - d^2 a_x^2 - d^2 a_y^2 + d^2 b_x^2 -3 d^2 b_y^2 + d^4\\ = 0 $