Subject: Problems (different) with solve and %solve
From: Emmanuel Charpentier
Date: Sun, 23 Dec 2012 16:05:20 +0000 (UTC)
Disclaimer : my trig is not merely rusty, it's corroded thru 'n thru (I
took it about 40 years ago and didn't use it much since : my main tool
for dealing with a root is a syndesmotome :-).
Maxima 5.28.0 http://maxima.sourceforge.net
using Lisp GNU Common Lisp (GCL) GCL 2.6.7 (a.k.a. GCL)
Distributed under the GNU Public License. See the file COPYING.
Dedicated to the memory of William Schelter.
The function bug_report() provides bug reporting information.
(%i1) load(to_poly_solve);
Loading maxima-grobner $Revision: 1.6 $ $Date: 2009-06-02 07:49:49 $
(%o1) /usr/share/maxima/5.28.0/share/to_poly_solve/to_poly_solve.mac
I saw this one in a Sage ticket about symbolic expression. It was some
years old, nd I wondered if this was still valid, t least as far Maxima
is concerned.
(%i2) foo:(sin(x) - 8*cos(x)*sin(x))*(sin(x)^2 + cos(x)) - (2*cos(x)*sin
(x) -sin(x))*(-2*sin(x)^2 + 2*cos(x)^2 - cos(x));
2
(%o2) (sin(x) - 8 cos(x) sin(x)) (sin (x) + cos(x))
2 2
- (2 cos(x) sin(x) - sin(x)) (- 2 sin (x) + 2 cos (x) -
cos(x))
(%i3) s1:solve(foo=0,x);
solve: using arc-trig functions to get a solution.
Some solutions will be lost.
%pi
(%o3) [x = %pi, x = ---, x = 0]
2
This solution appears fishy, and is easily contradicted by a graphical
representation (not shown here). Let's check it :
(%i4) map(lambda([s],subst(rhs(s),x,foo)),s1);
(%o4) [0, - 1, 0]
The second solution is false. What does to_poly_solve ?
(%i5) s2:nicedummies(%solve(foo=0,x));
3/2
(%o5) %union([x = 2 %pi %z0], [x = 2 %pi %z1 - atan(2 ) + %pi],
3/2
[x = 2 %pi %z2 + atan(2 ) - %pi], [x = 2 %pi %z3 +
%pi])
Let's check :
(%i6) z0:0, z1:0, z2:0, z3:0;
(%o6) 0
(%i7) makelist(subst(part(s2,i,1),foo),i,1,length(s2));
9/2 3/2 5/2 3/2
2 2 2 2
5 (---- + ----) 11 (- ---- - ----)
9 3 9 3
(%o7) [0, --------------- + ------------------,
9 9
9/2 3/2 5/2 3/2
2 2 2 2
5 (- ---- - ----) 11 (---- + ----
)
9 3 9 3
----------------- +
----------------, 0]
9 9
The numerical expression for the two last roots look strange to my
(untrained) eye, and are not *obviously* 0. But :
(%i8) %,numer;
(%o8) [0, 6.6613381477509392E-16, - 6.6613381477509392E-16, 0]
Seems reasonable enough.
Let's hunt the eoorr of solve() :
(%i9) e1:expand(foo);
3 3 3 2
(%o9) - 4 cos(x) sin (x) - sin (x) - 4 cos (x) sin(x) - 4 cos (x) sin(x)
(%i10) solve(e1=0,x);
solve: using arc-trig functions to get a solution.
Some solutions will be lost.
%pi
(%o10) [x = %pi, x = ---, x = 0]
2
Bummer...
(%i11) e2:trigexpand(e1);
3 3 3 2
(%o11) - 4 cos(x) sin (x) - sin (x) - 4 cos (x) sin(x) - 4 cos (x) sin(x)
(%i12) solve(e2=0,x);
solve: using arc-trig functions to get a solution.
Some solutions will be lost.
%pi
(%o12) [x = %pi, x = ---, x = 0]
2
Bummer again ...
(%i13) e3:factor(e2);
2 2 3 2
(%o13) - sin(x) (4 cos(x) sin (x) + sin (x) + 4 cos (x) + 4 cos (x))
(%i14) solve(e3,x);
solve: using arc-trig functions to get a solution.
Some solutions will be lost.
%pi
(%o14) [x = %pi, x = ---, x = 0]
2
Still bummer...
(%i15) e4:trigsimp(e3);
2
(%o15) (- 3 cos (x) - 4 cos(x) - 1) sin(x)
(%i16) s3:solve(e4=0,x);
solve: using arc-trig functions to get a solution.
Some solutions will be lost.
1
(%o16) [x = 0, x = %pi - acos(-), x = %pi]
3
Aha ! This looks interesting, but, prima facie, different from %solve()'s
solution. Let's check this :
(%i17) map(lambda([s],subst(rhs(s),x,foo)),s3);
9/2 3/2 5/2 3/2
2 2 2 2
5 (---- + ----) 11 (- ---- - ----)
9 3 9 3
(%o17) [0, --------------- + ------------------, 0]
9 9
The numerical expression for the last root is identical to the expression
of the third solution given by %solve() : co?ncidence ?
Numerical check :
(%i18) %,numer;
(%o18) [0, 6.6613381477509392E-16, 0]
Close enough. Let's explorethe two different expressin of the non-trivial
root(s). Numerically :
(%i19) acos(1/3)-atan(2^(3/2));
1 3/2
(%o19) acos(-) - atan(2 )
3
(%i20) %,numer;
(%o20) 0.0
Symbolically :
(%i21) is(equal(acos(1/3),atan(2^(3/2))));
(%o21) true
First conclusion : this demonstrates that my knowledge of trig is
insufficient to understand Maxima's problem with this expression :-).
By the way :
(%i22) trigsimp(foo);
2
(%o22) (- 3 cos (x) - 4 cos(x) - 1) sin(x)
(%i23) solve(%=0,x);
solve: using arc-trig functions to get a solution.
Some solutions will be lost.
1
(%o23) [x = 0, x = %pi - acos(-), x = %pi]
3
We could have avoided the e1, e2, e3, e4 expressions d go ttraigght
simplification.
Second conclusion : current Maxim's solve() has a bug that can be
triggered by some trigonometric expressions. This qualifies has a bug
because solve() gives a *false* answer without relevant wrning.
Third conclusion : %solve() gives a (seemingly) correct answer under a
complicated form, that Mxima can recognize as equivalent to the simpler
form given by solve() (when used on a simplified expression), but can't
regenerate itself :
(%i24) trigsimp(atan(2^(3/2)));
3/2
(%o24) atan(2 )
This qualifies as a wart, rather than a bug (the answer is, after all,
(seemingly) correct). It might derail further computations. However,
since there is no algorithmic definition of "simpler", I do not see a
soution to this wart.
HTH,
Emmanuel Charpentier