turning differences into sums

Hi all,

I have something like

(sqrt(a) - sqrt(b)) / c

and I would like to multiply to and bottom by

(sqrt(a) + sqrt(b))

and simplify each to get

(a - b) / (c * (sqrt(a) + sqrt(b)))

of course my actual expressions are much more complicated (and I hope my (a-b) will simplify further) and I cannot see how to do it programmatically without lots of manual copy and paste with potential for errors.

actual code:

s(x,y,z,w) := (x^2 + y^2 - z^2 - w^2) * (x^2 - y^2) / (x^2 + y^2);
t(x,y,z,w) := (x^2 + y^2 - z^2 - w^2) * (2 * x * y) / (x^2 + y^2);
u(x,y,z,w) := 2 * sqrt((x^2 + y^2) * (z^2 + w^2)) * (z^2 - w^2) / (z^2 + w^2);
v(x,y,z,w) := 2 * sqrt((x^2 + y^2) * (z^2 + w^2)) * (2 * z * w) / (z^2 + w^2);
p(f) := factor(expand(f(X+x,Y+y,Z+z,W+w) - f(X,Y,Z,W)));
p(s);
p(t);
p(u);
p(v);

the last two expressions are the ones I want to manipulate as described, here is the last one:

-((4*(W*sqrt(Y^2+X^2)*Z*sqrt(z^2+2*Z*z+w^2+2*W*w+Z^2+W^2)
            -sqrt(Z^2+W^2)*w*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)*z
            -W*sqrt(Z^2+W^2)*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)*z
            -Z*sqrt(Z^2+W^2)*w*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)
            -W*Z*sqrt(Z^2+W^2)*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)))
 /(sqrt(Z^2+W^2)*sqrt(z^2+2*Z*z+w^2+2*W*w+Z^2+W^2)))


context: I'm trying to apply perturbation techniques to iterations of the 4D Hopfbrot fractal, a cousin of the famous 3D Mandelbulb.  Perturbation techniques popularized in the last years allow computationally efficient deep zooming of 2D fractals like the Mandelbrot set, and recently have been applied to the Mandelbulb.  The main idea is to use "small" differences from one "large" reference orbit, and the key step is symbolically simplifying so there is no catastrophic cancellation when evaluating numerically (i.e. naively (X+x)-X will give 0 instead of x when x << X).

Thanks,


Claude

https://mathr.co.uk