Re: turning differences into sums

Stavros Macrakis <[email protected]> Sun, 21 Jun 2026 10:56:34 -0400
Newsgroups gmane.comp.mathematics.maxima.general
Message-ID <CACLVabXFyzgrCrO22BZ+pTRjZD5ftTJyVmGmj8d+1553izDqCg@mail.gmail.com>
--===============2615045256494987324==
Content-Type: multipart/alternative; boundary="0000000000002ad8b00654c4bce4"

--0000000000002ad8b00654c4bce4
Content-Type: text/plain; charset="UTF-8"
Content-Transfer-Encoding: quoted-printable

Let's start with p(v):

 -((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)))

and extract the relevant part:

p111: part(p(v),1,1,1)

By inspection, the positive part is *first(p111)* and the negative part is
*rest(p111)*.

But to be more systematic, you could do something like this:

syntactically_negative_p(ex):=3D
    if numberp(ex)
      then is(ex<0)
    else not mapatom(ex) and is(op(ex)=3D"-");

partition_sum(expr,pred):=3D
  block([op,args,yes:[],no:[]],
    if mapatom(expr) or op(expr)#"+"
      then (op:"+",args:[expr])
    else (op:op(expr),args:args(expr)),
    for arg in args do
       if apply(pred,[arg])
         then push(arg,yes)
       else push(arg,no),
    [funmake(op,reverse(yes)),
     funmake(op,reverse(no))]
    )$

partition_sum(p111,syntactically_negative_p) =3D>
    [-(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),
           W*sqrt(Y^2+X^2)*Z*sqrt(z^2+2*Z*z+w^2+2*W*w+Z^2+W^2)]




On Sun, Jun 21, 2026 at 10:01=E2=80=AFAM Claude Heiland-Allen <claude@mathr=
.co.uk>
wrote:

> perfect! Thanks Stavros.
>
>
> here's that trick applied to my problem:
>
>
> display2d : false;
>
> s(x,y,z,w) :=3D (x^2 + y^2 - z^2 - w^2) * (x^2 - y^2) / (x^2 + y^2);
>
> t(x,y,z,w) :=3D (x^2 + y^2 - z^2 - w^2) * (2 * x * y) / (x^2 + y^2);
>
> u(x,y,z,w) :=3D 2 * sqrt((x^2 + y^2) * (z^2 + w^2)) * (z^2 - w^2) / (z^2 =
+
> w^2);
>
> v(x,y,z,w) :=3D 2 * sqrt((x^2 + y^2) * (z^2 + w^2)) * (2 * z * w) / (z^2 =
+
> w^2);
>
> p(f) :=3D factor(expand(f(X+x,Y+y,Z+z,W+w) - f(X,Y,Z,W)));
>                                                                       p(s=
);
> p(t);
> p(u);
>  p(v);
>
>                                                    m : substpart("-",
> expand(args(-p(v))[1]), 0);
>
> apply("/", [-expand(args(-p(v))[1] * m), ratsimp(args(-p(v))[2] * m)]);
>
>
>
>
> Thanks,
>
>
> Claude
>
>
> On 21 June 2026 14:24:54 CEST, Stavros Macrakis <[email protected]>
> wrote:
> >BTW, a trick for getting A+B from A-B is
> >
> >substpart("-",ex,0)
> >
> >Why - and not +? Because Maxima treats A-B as A+(-B)
> >
> >
> >On Sun, Jun 21, 2026, 07:39 Barton Willis <[email protected]> wrote:
> >
> >> If your primary goal is accurate floating-point evaluation, you might
> try
> >> Herbie <https://github.com/herbie-fp/herbie>.  But the online
> >> <http://herbie.uwplse.org> interface to Herbie timed out when I tried
> >> loading it a moment ago.
> >>
> >> Let us know what you all discover.
> >>
> >> --Barton
> >>
> >>
> >>
> >>
> >>
> >>
> >> ------------------------------
> >> *From:* Claude Heiland-Allen <[email protected]>
> >> *Sent:* Sunday, June 21, 2026 6:02 AM
> >> *To:* Stavros Macrakis <[email protected]>
> >> *Cc:* maxima-discuss <[email protected]>
> >> *Subject:* Re: [Maxima-discuss] turning differences into sums
> >>
> >> Caution: Non-NU Email
> >>
> >>
> >> Thanks for your reply Stavros, however the first example doesn't make
> >> any difference in this case and the second one has the copy/paste issu=
e
> >> that I'm trying to avoid (though, it uses a nicer syntax than my
> >> previous attempts, which is good to learn).
> >>
> >> Concretely, I want to programmatically transform all the subtractions =
in
> >> the numerator to additions, that is, go from
> >>
> >> -((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)))
> >>
> >>
> >> to
> >>
> >>
> >> -((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)))
> >>
> >>
> >> so that I can the multiply the parts of the fraction by it.
> >>
> >>
> >> Thanks,
> >>
> >>
> >> Claude
> >>
> >> On 21/06/2026 12:49, Stavros Macrakis wrote:
> >> > Try
> >> >
> >> > factor(1/ratsimp(1/ex)),algebraic;
> >> >
> >> > ratsimp/algebraic generally eliminates roots in the denominator
> >> >
> >> > or more "manually"
> >> >
> >> > apply("/",factor(expand(args(ex)*(sqrt(a)+sqrt(b)))))
> >> >
> >> >
> >> > On Sun, Jun 21, 2026, 04:46 Claude Heiland-Allen <[email protected]=
>
> >> > wrote:
> >> >
> >> >     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) :=3D (x^2 + y^2 - z^2 - w^2) * (x^2 - y^2) / (x^2 + y=
^2);
> >> >     t(x,y,z,w) :=3D (x^2 + y^2 - z^2 - w^2) * (2 * x * y) / (x^2 + y=
^2);
> >> >     u(x,y,z,w) :=3D 2 * sqrt((x^2 + y^2) * (z^2 + w^2)) * (z^2 - w^2=
) /
> >> >     (z^2 + w^2);
> >> >     v(x,y,z,w) :=3D 2 * sqrt((x^2 + y^2) * (z^2 + w^2)) * (2 * z * w=
) /
> >> >     (z^2 + w^2);
> >> >     p(f) :=3D 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 iteratio=
ns
> >> >     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 numericall=
y
> >> >     (i.e. naively (X+x)-X will give 0 instead of x when x << X).
> >> >
> >> >     Thanks,
> >> >
> >> >
> >> >     Claude
> >> >
> >> >
> >>
> https://urldefense.com/v3/__https://mathr.co.uk__;!!PvXuogZ4sRB2p-tU!HNBM=
Y8w_rKA70J6zDmTJj3M67tGC5xyBsRZ4uBzQTT764MBPm7zFlVo4-1aaHkqqukTjwTljmsTidt0=
$
> >> >
> >> >
> >> >     _______________________________________________
> >> >     Maxima-discuss mailing list
> >> >     [email protected]
> >> >
> >>
> https://urldefense.com/v3/__https://lists.sourceforge.net/lists/listinfo/=
maxima-discuss__;!!PvXuogZ4sRB2p-tU!HNBMY8w_rKA70J6zDmTJj3M67tGC5xyBsRZ4uBz=
QTT764MBPm7zFlVo4-1aaHkqqukTjwTljtm35Guw$
> >> >
> >>
> >>
> >> _______________________________________________
> >> Maxima-discuss mailing list
> >> [email protected]
> >>
> >>
> https://urldefense.com/v3/__https://lists.sourceforge.net/lists/listinfo/=
maxima-discuss__;!!PvXuogZ4sRB2p-tU!HNBMY8w_rKA70J6zDmTJj3M67tGC5xyBsRZ4uBz=
QTT764MBPm7zFlVo4-1aaHkqqukTjwTljtm35Guw$
> >>
>
> https://mathr.co.uk
>

--0000000000002ad8b00654c4bce4
Content-Type: text/html; charset="UTF-8"
Content-Transfer-Encoding: quoted-printable

<div dir=3D"ltr"><div class=3D"gmail_default" style=3D"font-family:georgia,=
serif;font-size:small;color:#330000">Let&#39;s start with p(v):<br><br></di=
v><div class=3D"gmail_default" style=3D"font-size:small;color:rgb(51,0,0)">=
<span style=3D"font-family:georgia,serif">=C2=A0</span><font face=3D"monosp=
ace">-((4*(W*sqrt(Y^2+X^2)*Z*sqrt(z^2+2*Z*z+w^2+2*W*w+Z^2+W^2)<br>=C2=A0 =
=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 -sqrt(Z^2+W^2)*w*sqrt(y^2+2*Y*y+x^2+2*X*=
x+Y^2+X^2)*z<br>=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 -W*sqrt(Z^2+W^2)*=
sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)*z<br>=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =
=C2=A0 -Z*sqrt(Z^2+W^2)*w*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)<br>=C2=A0 =C2=
=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 -W*Z*sqrt(Z^2+W^2)*sqrt(y^2+2*Y*y+x^2+2*X*x=
+Y^2+X^2)))<br>=C2=A0/(sqrt(Z^2+W^2)*sqrt(z^2+2*Z*z+w^2+2*W*w+Z^2+W^2)))</f=
ont></div><div class=3D"gmail_default" style=3D"font-family:georgia,serif;f=
ont-size:small;color:#330000"><br></div><div class=3D"gmail_default" style=
=3D"font-family:georgia,serif;font-size:small;color:#330000">and extract th=
e relevant part:</div><div class=3D"gmail_default" style=3D"font-family:geo=
rgia,serif;font-size:small;color:#330000"><br></div><div class=3D"gmail_def=
ault" style=3D"font-size:small;color:rgb(51,0,0)"><font face=3D"monospace">=
p111: part(p(v),1,1,1)</font></div><div class=3D"gmail_default" style=3D"fo=
nt-family:georgia,serif;font-size:small;color:#330000"><br></div><div class=
=3D"gmail_default" style=3D"font-family:georgia,serif;font-size:small;color=
:#330000">By inspection, the positive part is <b>first(p111)</b> and the ne=
gative part is <b>rest(p111)</b>.</div><div class=3D"gmail_default" style=
=3D"font-family:georgia,serif;font-size:small;color:#330000"><br></div><div=
 class=3D"gmail_default" style=3D"font-family:georgia,serif;font-size:small=
;color:#330000">But to be more systematic, you could do something like this=
:</div><div class=3D"gmail_default" style=3D"font-family:georgia,serif;font=
-size:small;color:#330000"><br></div><div class=3D"gmail_default" style=3D"=
font-size:small;color:rgb(51,0,0)"><font face=3D"monospace">syntactically_n=
egative_p(ex):=3D<br>=C2=A0 =C2=A0 if numberp(ex)<br>=C2=A0 =C2=A0 =C2=A0 t=
hen is(ex&lt;0)<br>=C2=A0 =C2=A0 else not mapatom(ex) and is(op(ex)=3D&quot=
;-&quot;);<br><br>partition_sum(expr,pred):=3D<br>=C2=A0 block([op,args,yes=
:[],no:[]],<br>=C2=A0 =C2=A0 if mapatom(expr) or op(expr)#&quot;+&quot;<br>=
=C2=A0 =C2=A0 =C2=A0 then (op:&quot;+&quot;,args:[expr])<br>=C2=A0 =C2=A0 e=
lse (op:op(expr),args:args(expr)),<br>=C2=A0 =C2=A0 for arg in args do<br>=
=C2=A0 =C2=A0 =C2=A0 =C2=A0if apply(pred,[arg])<br>=C2=A0 =C2=A0 =C2=A0 =C2=
=A0 =C2=A0then push(arg,yes)<br>=C2=A0 =C2=A0 =C2=A0 =C2=A0else push(arg,no=
),<br>=C2=A0 =C2=A0 [funmake(op,reverse(yes)),<br>=C2=A0 =C2=A0 =C2=A0funma=
ke(op,reverse(no))]<br>=C2=A0 =C2=A0 )$<br></font></div><div class=3D"gmail=
_default" style=3D"font-size:small;color:rgb(51,0,0)"><font face=3D"monospa=
ce"><br></font></div><div class=3D"gmail_default" style=3D"font-size:small;=
color:rgb(51,0,0)"><font face=3D"monospace">partition_sum(p111,syntacticall=
y_negative_p) =3D&gt;</font></div><div class=3D"gmail_default" style=3D"fon=
t-size:small;color:rgb(51,0,0)"><font face=3D"monospace">=C2=A0 =C2=A0<span=
 style=3D"background-color:transparent">=C2=A0</span><span style=3D"backgro=
und-color:transparent;color:rgb(34,34,34)">[-(sqrt(Z^2+W^2)*w*sqrt(y^2+2*Y*=
y+x^2+2*X*x+Y^2+X^2)*z)</span></font></div><div class=3D"gmail_default" sty=
le=3D"font-size:small;color:rgb(51,0,0)"><font face=3D"monospace"><span sty=
le=3D"background-color:transparent;color:rgb(34,34,34)">=C2=A0 =C2=A0 =C2=
=A0-W*sqrt(Z^2+W^2)*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)*z</span></font></div>=
<div class=3D"gmail_default" style=3D"font-size:small;color:rgb(51,0,0)"><f=
ont face=3D"monospace"><span style=3D"background-color:transparent;color:rg=
b(34,34,34)">=C2=A0 =C2=A0 =C2=A0-Z*sqrt(Z^2+W^2)*w*sqrt(y^2+2*Y*y+x^2+2*X*=
x+Y^2+X^2)</span></font></div><font face=3D"monospace">=C2=A0 =C2=A0 =C2=A0=
-W*Z*sqrt(Z^2+W^2)*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2),</font><div><font face=
=3D"monospace"><span class=3D"gmail_default" style=3D"font-size:small;color=
:rgb(51,0,0)">=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0</span>W*sqrt(Y^2+X^=
2)*Z*sqrt(z^2+2*Z*z+w^2+2*W*w+Z^2+W^2)]</font></div><div><font face=3D"mono=
space"><br></font><div class=3D"gmail_default" style=3D"font-size:small;col=
or:rgb(51,0,0)"><br></div><div class=3D"gmail_default" style=3D"font-family=
:georgia,serif;font-size:small;color:#330000"><br></div></div></div><br><di=
v class=3D"gmail_quote gmail_quote_container"><div dir=3D"ltr" class=3D"gma=
il_attr">On Sun, Jun 21, 2026 at 10:01=E2=80=AFAM Claude Heiland-Allen &lt;=
<a href=3D"mailto:[email protected]">[email protected]</a>&gt; wrote:<br>=
</div><blockquote class=3D"gmail_quote" style=3D"margin:0px 0px 0px 0.8ex;b=
order-left:1px solid rgb(204,204,204);padding-left:1ex">perfect! Thanks Sta=
vros.<br>
<br>
<br>
here&#39;s that trick applied to my problem:<br>
<br>
<br>
display2d : false;=C2=A0 =C2=A0 =C2=A0 =C2=A0<br>
<br>
s(x,y,z,w) :=3D (x^2 + y^2 - z^2 - w^2) * (x^2 - y^2) / (x^2 + y^2);=C2=A0 =
=C2=A0 <br>
<br>
t(x,y,z,w) :=3D (x^2 + y^2 - z^2 - w^2) * (2 * x * y) / (x^2 + y^2);<br>
<br>
u(x,y,z,w) :=3D 2 * sqrt((x^2 + y^2) * (z^2 + w^2)) * (z^2 - w^2) / (z^2 + =
w^2);<br>
<br>
v(x,y,z,w) :=3D 2 * sqrt((x^2 + y^2) * (z^2 + w^2)) * (2 * z * w) / (z^2 + =
w^2);=C2=A0 =C2=A0 =C2=A0 =C2=A0 <br>
<br>
p(f) :=3D factor(expand(f(X+x,Y+y,Z+z,W+w) - f(X,Y,Z,W)));<br>
=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=
=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =
=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=
=A0 =C2=A0 =C2=A0 =C2=A0 p(s);<br>
p(t);<br>
p(u);=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =
=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=
=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =
=C2=A0 =C2=A0p(v);=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 <br>
=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=
=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =
=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=
=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =
=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=
=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0m =
: substpart(&quot;-&quot;, expand(args(-p(v))[1]), 0);<br>
<br>
apply(&quot;/&quot;, [-expand(args(-p(v))[1] * m), ratsimp(args(-p(v))[2] *=
 m)]);=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0=
 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=
=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =
=C2=A0 =C2=A0 =C2=A0 =C2=A0<br>
<br>
<br>
<br>
Thanks,<br>
<br>
<br>
Claude<br>
<br>
<br>
On 21 June 2026 14:24:54 CEST, Stavros Macrakis &lt;<a href=3D"mailto:macra=
[email protected]" target=3D"_blank">[email protected]</a>&gt; wrote:<br>
&gt;BTW, a trick for getting A+B from A-B is<br>
&gt;<br>
&gt;substpart(&quot;-&quot;,ex,0)<br>
&gt;<br>
&gt;Why - and not +? Because Maxima treats A-B as A+(-B)<br>
&gt;<br>
&gt;<br>
&gt;On Sun, Jun 21, 2026, 07:39 Barton Willis &lt;<a href=3D"mailto:willisb=
@unk.edu" target=3D"_blank">[email protected]</a>&gt; wrote:<br>
&gt;<br>
&gt;&gt; If your primary goal is accurate floating-point evaluation, you mi=
ght try<br>
&gt;&gt; Herbie &lt;<a href=3D"https://github.com/herbie-fp/herbie" rel=3D"=
noreferrer" target=3D"_blank">https://github.com/herbie-fp/herbie</a>&gt;.=
=C2=A0 But the online<br>
&gt;&gt; &lt;<a href=3D"http://herbie.uwplse.org" rel=3D"noreferrer" target=
=3D"_blank">http://herbie.uwplse.org</a>&gt; interface to Herbie timed out =
when I tried<br>
&gt;&gt; loading it a moment ago.<br>
&gt;&gt;<br>
&gt;&gt; Let us know what you all discover.<br>
&gt;&gt;<br>
&gt;&gt; --Barton<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; ------------------------------<br>
&gt;&gt; *From:* Claude Heiland-Allen &lt;<a href=3D"mailto:[email protected]=
.uk" target=3D"_blank">[email protected]</a>&gt;<br>
&gt;&gt; *Sent:* Sunday, June 21, 2026 6:02 AM<br>
&gt;&gt; *To:* Stavros Macrakis &lt;<a href=3D"mailto:[email protected]" t=
arget=3D"_blank">[email protected]</a>&gt;<br>
&gt;&gt; *Cc:* maxima-discuss &lt;<a href=3D"mailto:[email protected]=
urceforge.net" target=3D"_blank">[email protected]</a>&g=
t;<br>
&gt;&gt; *Subject:* Re: [Maxima-discuss] turning differences into sums<br>
&gt;&gt;<br>
&gt;&gt; Caution: Non-NU Email<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; Thanks for your reply Stavros, however the first example doesn&#39=
;t make<br>
&gt;&gt; any difference in this case and the second one has the copy/paste =
issue<br>
&gt;&gt; that I&#39;m trying to avoid (though, it uses a nicer syntax than =
my<br>
&gt;&gt; previous attempts, which is good to learn).<br>
&gt;&gt;<br>
&gt;&gt; Concretely, I want to programmatically transform all the subtracti=
ons in<br>
&gt;&gt; the numerator to additions, that is, go from<br>
&gt;&gt;<br>
&gt;&gt; -((4*(W*sqrt(Y^2+X^2)*Z*sqrt(z^2+2*Z*z+w^2+2*W*w+Z^2+W^2)<br>
&gt;&gt;=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 -sqrt(Z^2+W^2)*w*s=
qrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)*z<br>
&gt;&gt;=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 -W*sqrt(Z^2+W^2)*s=
qrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)*z<br>
&gt;&gt;=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0 -Z*sqrt(Z^2+W^2)*w=
*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)<br>
&gt;&gt;<br>
&gt;&gt; -W*Z*sqrt(Z^2+W^2)*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)))<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; to<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; -((4*(W*sqrt(Y^2+X^2)*Z*sqrt(z^2+2*Z*z+w^2+2*W*w+Z^2+W^2)<br>
&gt;&gt; +sqrt(Z^2+W^2)*w*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)*z<br>
&gt;&gt;<br>
&gt;&gt; +W*sqrt(Z^2+W^2)*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)*z<br>
&gt;&gt; +Z*sqrt(Z^2+W^2)*w*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)<br>
&gt;&gt; +W*Z*sqrt(Z^2+W^2)*sqrt(y^2+2*Y*y+x^2+2*X*x+Y^2+X^2)))<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; so that I can the multiply the parts of the fraction by it.<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; Thanks,<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; Claude<br>
&gt;&gt;<br>
&gt;&gt; On 21/06/2026 12:49, Stavros Macrakis wrote:<br>
&gt;&gt; &gt; Try<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt; factor(1/ratsimp(1/ex)),algebraic;<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt; ratsimp/algebraic generally eliminates roots in the denominat=
or<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt; or more &quot;manually&quot;<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt; apply(&quot;/&quot;,factor(expand(args(ex)*(sqrt(a)+sqrt(b)))=
))<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt; On Sun, Jun 21, 2026, 04:46 Claude Heiland-Allen &lt;<a href=
=3D"mailto:[email protected]" target=3D"_blank">[email protected]</a>&gt;=
<br>
&gt;&gt; &gt; wrote:<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0Hi all,<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0I have something like<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0(sqrt(a) - sqrt(b)) / c<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0and I would like to multiply to and bottom=
 by<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0(sqrt(a) + sqrt(b))<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0and simplify each to get<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0(a - b) / (c * (sqrt(a) + sqrt(b)))<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0of course my actual expressions are much m=
ore complicated (and I<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0hope my (a-b) will simplify further) and I=
 cannot see how to do it<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0programmatically without lots of manual co=
py and paste with<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0potential for errors.<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0actual code:<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0s(x,y,z,w) :=3D (x^2 + y^2 - z^2 - w^2) * =
(x^2 - y^2) / (x^2 + y^2);<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0t(x,y,z,w) :=3D (x^2 + y^2 - z^2 - w^2) * =
(2 * x * y) / (x^2 + y^2);<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0u(x,y,z,w) :=3D 2 * sqrt((x^2 + y^2) * (z^=
2 + w^2)) * (z^2 - w^2) /<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0(z^2 + w^2);<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0v(x,y,z,w) :=3D 2 * sqrt((x^2 + y^2) * (z^=
2 + w^2)) * (2 * z * w) /<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0(z^2 + w^2);<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0p(f) :=3D factor(expand(f(X+x,Y+y,Z+z,W+w)=
 - f(X,Y,Z,W)));<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0p(s);<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0p(t);<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0p(u);<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0p(v);<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0the last two expressions are the ones I wa=
nt to manipulate as<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0described, here is the last one:<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0-((4*(W*sqrt(Y^2+X^2)*Z*sqrt(z^2+2*Z*z+w^2=
+2*W*w+Z^2+W^2)<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0-sqrt(Z^2+W^2)*w*sqrt(y^2+2*Y*y+x^2+2*X*x+=
Y^2+X^2)*z<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0-W*sqrt(Z^2+W^2)*sqrt(y^2+2*Y*y+x^2+2*X*x+=
Y^2+X^2)*z<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0-Z*sqrt(Z^2+W^2)*w*sqrt(y^2+2*Y*y+x^2+2*X*=
x+Y^2+X^2)<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0-W*Z*sqrt(Z^2+W^2)*sqrt(y^2+2*Y*y+x^2+2*X*=
x+Y^2+X^2)))<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0 /(sqrt(Z^2+W^2)*sqrt(z^2+2*Z*z+w^2+2*W*w+=
Z^2+W^2)))<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0context: I&#39;m trying to apply perturbat=
ion techniques to iterations<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0of the 4D Hopfbrot fractal, a cousin of th=
e famous 3D Mandelbulb.<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0Perturbation techniques popularized in the=
 last years allow<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0computationally efficient deep zooming of =
2D fractals like the<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0Mandelbrot set, and recently have been app=
lied to the Mandelbulb.<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0The main idea is to use &quot;small&quot; =
differences from one &quot;large&quot;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0reference orbit, and the key step is symbo=
lically simplifying so<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0there is no catastrophic cancellation when=
 evaluating numerically<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0(i.e. naively (X+x)-X will give 0 instead =
of x when x &lt;&lt; X).<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0Thanks,<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0Claude<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;<br>
&gt;&gt; <a href=3D"https://urldefense.com/v3/__https://mathr.co.uk__;!!PvX=
uogZ4sRB2p-tU!HNBMY8w_rKA70J6zDmTJj3M67tGC5xyBsRZ4uBzQTT764MBPm7zFlVo4-1aaH=
kqqukTjwTljmsTidt0$" rel=3D"noreferrer" target=3D"_blank">https://urldefens=
e.com/v3/__https://mathr.co.uk__;!!PvXuogZ4sRB2p-tU!HNBMY8w_rKA70J6zDmTJj3M=
67tGC5xyBsRZ4uBzQTT764MBPm7zFlVo4-1aaHkqqukTjwTljmsTidt0$</a><br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0__________________________________________=
_____<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0Maxima-discuss mailing list<br>
&gt;&gt; &gt;=C2=A0 =C2=A0 =C2=A0<a href=3D"mailto:[email protected]=
rceforge.net" target=3D"_blank">[email protected]</a><br=
>
&gt;&gt; &gt;<br>
&gt;&gt; <a href=3D"https://urldefense.com/v3/__https://lists.sourceforge.n=
et/lists/listinfo/maxima-discuss__;!!PvXuogZ4sRB2p-tU!HNBMY8w_rKA70J6zDmTJj=
3M67tGC5xyBsRZ4uBzQTT764MBPm7zFlVo4-1aaHkqqukTjwTljtm35Guw$" rel=3D"norefer=
rer" target=3D"_blank">https://urldefense.com/v3/__https://lists.sourceforg=
e.net/lists/listinfo/maxima-discuss__;!!PvXuogZ4sRB2p-tU!HNBMY8w_rKA70J6zDm=
TJj3M67tGC5xyBsRZ4uBzQTT764MBPm7zFlVo4-1aaHkqqukTjwTljtm35Guw$</a><br>
&gt;&gt; &gt;<br>
&gt;&gt;<br>
&gt;&gt;<br>
&gt;&gt; _______________________________________________<br>
&gt;&gt; Maxima-discuss mailing list<br>
&gt;&gt; <a href=3D"mailto:[email protected]" target=3D"=
_blank">[email protected]</a><br>
&gt;&gt;<br>
&gt;&gt; <a href=3D"https://urldefense.com/v3/__https://lists.sourceforge.n=
et/lists/listinfo/maxima-discuss__;!!PvXuogZ4sRB2p-tU!HNBMY8w_rKA70J6zDmTJj=
3M67tGC5xyBsRZ4uBzQTT764MBPm7zFlVo4-1aaHkqqukTjwTljtm35Guw$" rel=3D"norefer=
rer" target=3D"_blank">https://urldefense.com/v3/__https://lists.sourceforg=
e.net/lists/listinfo/maxima-discuss__;!!PvXuogZ4sRB2p-tU!HNBMY8w_rKA70J6zDm=
TJj3M67tGC5xyBsRZ4uBzQTT764MBPm7zFlVo4-1aaHkqqukTjwTljtm35Guw$</a><br>
&gt;&gt;<br>
<br>
<a href=3D"https://mathr.co.uk" rel=3D"noreferrer" target=3D"_blank">https:=
//mathr.co.uk</a><br>
</blockquote></div>

--0000000000002ad8b00654c4bce4--


--===============2615045256494987324==
Content-Type: text/plain; charset="us-ascii"
MIME-Version: 1.0
Content-Transfer-Encoding: 7bit
Content-Disposition: inline


--===============2615045256494987324==
Content-Type: text/plain; charset="us-ascii"
MIME-Version: 1.0
Content-Transfer-Encoding: 7bit
Content-Disposition: inline

_______________________________________________
Maxima-discuss mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/maxima-discuss

--===============2615045256494987324==--