The option variable realonly is quite confusing since it does not find
real solutions, but only solutions which are free of %i (purely using
freeof).
For example, when realonly is false,
algsys ([x^4 + 1], [x]);
returns
[[x = (-1)^(1/4)], [x = -(-1)^(1/4)*%i],
[x = -(-1)^(1/4)], [x = (-1)^(1/4)*%i]]
But when realonly is true, it returns
[[x = (-1)^(1/4)], [x = -(-1)^(1/4)]]
while it is natural to expect [].
Maybe it is better to modify the behavior and filter roots by
checking something like
is_real(x) := is(trigsimp(imagpart(x)) = 0);
and not just
freeof(%i, x)
With the freeof approach we may also omit real roots of irreducible
polynomials:
sol : map ('rhs, solve (3*x^3 - 3*x + 1));
[(sqrt(3)*%i/2-1/2)/(3*(3^-(3/2)*%i/2-1/6)^(1/3))
+(3^-(3/2)*%i/2-1/6)^(1/3)*(-sqrt(3)*%i/2-1/2),
(3^-(3/2)*%i/2-1/6)^(1/3)*(sqrt(3)*%i/2-1/2)
+(-sqrt(3)*%i/2-1/2)/(3*(3^-(3/2)*%i/2-1/6)^(1/3)),
(3^-(3/2)*%i/2-1/6)^(1/3)+1/(3*(3^-(3/2)*%i/2-1/6)^(1/3))]
map (lambda([x], freeof(%i, x)), sol);
[false,false,false]
map ('is_real, sol);
[true,true,true]
So when realonly and algexact are set to true,
algsys ([3*x^3 - 3*x + 1], [x])
just returns [].
Volker van Nek
2010-01-23
Hi Alexey,
yesterday I noticed this problem too and found your bug report.
I think you are right, freeof(%i, x) is not a sufficient test to recognize non-zero imaginary parts.
E.g. (-1)^(1/4) passes.
So I wonder why you did not fix that problem. You suggested the right thing.
It is one line in src/algsys.lisp to change:
348-352:
(defun realonly (rootsl)
(cond ((null rootsl) nil)
;; ((freeof '$%i (car rootsl)) ;; problem
((equal 0 (sratsimp ($imagpart (car rootsl)))) ;; fix ?
(nconc (list (car rootsl)) (realonly (cdr rootsl))))
(t (realonly (cdr rootsl)))))
This change passes the test suite. So I would commit this new line.
Or do I overlook something? Why did you hesitate to commit?
Volker van Nek
Volker van Nek
2010-03-27
my recent cvs commit closes this report
Volker van Nek
Volker van Nek
2010-03-27