[Maxima-bugs] [ maxima-Bugs-2184396 ] Wrong factorization of sqrt()

SourceForge Headquarters 225 Broadway Suite 1600 San Diego, CA 92101 +1 (858) 454-5900

Bugs item #2184396, was opened at 2008-10-21 15:05
Message generated for change (Comment added) made by crategus
You can respond by visiting: 
https://sourceforge.net/tracker/?func=detail&atid=104933&aid=2184396&group_id=4933

Please note that this message will contain a full copy of the comment thread,
including the initial issue submission, for this request,
not just the latest update.
Category: None
Group: None
Status: Open
Resolution: None
Priority: 5
Private: No
Submitted By: Satoshi Adachi (satoshi_adachi)
Assigned to: Nobody/Anonymous (nobody)
Summary: Wrong factorization of sqrt()

Initial Comment:
Dear Developers of Maxima,

I found a wrong behavihor of sqrt().
It is a principle that sqrt(X*Y) should not be factorized to sqrt(X)*sqrt(Y)
unless X and Y can be judged to be non-negative.
I found a case in which this principle is broken.

Here is a program for demonstration:
------------------------------------------------------------------------------
/*
 *   wrong_factorization_of_sqrt.maxima:
 *
 *   S.Adachi 2008/10/01
 */

display2d:false;

/* real for 1 <= x <= 2 */
correct:sqrt((1-1/x)*(2/x-1));

/* real for 2-sqrt(2) <= x <= 2+sqrt(2) */
INCORRECT:sqrt((1-(2-sqrt(2))/x)*((2+sqrt(2))/x-1));

INCORRECT_AT_2:at(INCORRECT,[x=2]);

SHOULD_BE_POSITIVE:float(INCORRECT_AT_2);

/* END */
------------------------------------------------------------------------------

The result of execution is as follows:
------------------------------------------------------------------------------
Maxima 5.14.0cvs http://maxima.sourceforge.net
Using Lisp GNU Common Lisp (GCL) GCL 2.6.7 (aka 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)              batch(wrong_factorization_of_sqrt.maxima)

batching #p/Volumes/HFS+2/home/adachi/work/307/wrong_factorization_of_sqrt.maxima
(%i2)                          display2d : false
(%o2) false
(%i3) correct:sqrt((1-1/x)*(2/x-1))
(%o3) sqrt((1-1/x)*(2/x-1))
(%i4) INCORRECT:sqrt((1-(2-sqrt(2))/x)*((sqrt(2)+2)/x-1))
(%o4) sqrt((2-sqrt(2))/x-1)*sqrt(1-(sqrt(2)+2)/x)
(%i5) INCORRECT_AT_2:at(INCORRECT,[x = 2])
(%o5) sqrt((2-sqrt(2))/2-1)*sqrt(1-(sqrt(2)+2)/2)
(%i6) SHOULD_BE_POSITIVE:float(INCORRECT_AT_2)
(%o6) -0.70710678118655
(%o7) "wrong_factorization_of_sqrt.maxima"
------------------------------------------------------------------------------

At first, as is seen from (%i2) and (%o2),
sqrt(1-1/x)*(2/x-1)) is not factorized to sqrt(...)*sqrt(...)
in any automatic way. This is the correct behavihor.

Next, as is seen from (%i3) and (%o3),
sqrt(1-(2-sqrt(2))/x)*((2+sqrt(2))/x-1)) is factorized to
sqrt((2-sqrt(2))/x-1)*sqrt(1-(sqrt(2)+2)/x) in an automatic way.
THIS IS AN INCORRECT BEHAVIOR.

This expression is non-negative for 2-sqrt(2) <= x <= 2+sqrt(2).
However, as is demonstrated by (%i4)-(%o5),
the calculated expression takes a negative value at x=2.

Please stop the incorrect factorization of sqrt(X*Y) to sqrt(X)*sqrt(Y).

Sincerely yours,

Satoshi Adachi

----------------------------------------------------------------------

>Comment By: Dieter Kaiser (crategus)
Date: 2009-11-11 23:43

Message:
I think I have found the bug in the routine compsplt.

(defun compsplt (x)
  (cond ((atom x) (setq lhs x rhs 0))
	((atom (car x)) (setq lhs x rhs 0))
--->	((not (null (cdr (symbols x)))) (compsplt2 x))
	(t (compsplt1 x))))

At the marked line Maxima tests the CDR of the result of the function
SYMBOLS. SYMBOLS returns a list of all variables of an expression. The list
is of the form '($X $Y ...), where $X, $Y are the variables. Because
COMPSPLT takes the CDR of SYMBOLS, one symbol in an expression is not
recognized. Later in the algorithm the sign of this variable is not taken
into account.

For almost all expression this behavior seems to be no problem, but not
for the examples of this bug report.

When we change the test of the marked line:

   ((not (null (symbols x))) ...

the sign of the given examples of this bug report are correct:

(%i3) sign(1-(2+sqrt(2))*x);
(%o3) pnz
(%i4) sign(1-(2+sqrt(2))/x);
(%o4) pnz

and the reported expression does not longer factorize wrongly:

(%i5) sqrt((1-(2-sqrt(2))/x)*((2+sqrt(2))/x-1));
(%o5) sqrt((1-(2-sqrt(2))/x)*((sqrt(2)+2)/x-1))

One example of the testsuite no longer works. It is one of the new
examples for the sign of the abs function. This is the example:

(%i2) assume(abs(x)<%e/2+sin(1));
(%o2) [sin(1)+%e/2 > abs(x)]
(%i3)  is(x>-2*%e);
(%o3) true

This will be the new result:

(%i6) assume(abs(x)<%e/2+sin(1));
(%o6) [-abs(x)+sin(1)+%e/2 > 0]

(%i7) is(x>-2*%e);
(%o7) unknown

This test no longer works, because the ordering of the terms will change.
The reason is that COMPSPLT calls splitsum for one variable x in the
expression too.

So, with one exception the testsuite and the share_testsuite have no
problems with the suggested change in COMPSPLT.

Dieter Kaiser

----------------------------------------------------------------------

Comment By: Dieter Kaiser (crategus)
Date: 2009-11-11 01:27

Message:
For the record. Again two example for a wrong sign:

(%i3) sign(1-(1+sqrt(2))/x);
(%o3) neg
(%i4) sign(1-(1+sqrt(2))*x);
(%o4) neg

For these expressions the functions splitprod and splitsum are called. I
have not figured out up to now the reason for the bug, but the variable x
is lost somewhere and not taken into account to determine the sign.

Dieter Kaiser

----------------------------------------------------------------------

Comment By: Dieter Kaiser (crategus)
Date: 2009-10-08 00:10

Message:
Maxima can simplify the sqrt function correctly.

No simplification, when the sign of a is not known:

(%i1) sqrt(a*z);
(%o1) sqrt(a*z)

Simplification for a a positive or b a negative value:

(%i2) assume(a>0)$
(%i3) sqrt(a*z);
(%o3) sqrt(a)*sqrt(z)

(%i4) assume(b<0)$
(%i5) sqrt(b*z);
(%o5) sqrt(-b)*sqrt(-z)

The problem of this bug report is one of the factors of the example. The
sign is wrongly determined to be negative:

(%i8) sign(1-(2-sqrt(2))/x);
(%o8) neg

Therefore Maxima does the above simplification, which is wrong.

Dieter Kaiser

----------------------------------------------------------------------

Comment By: boud  (boud1)
Date: 2009-02-19 03:31

Message:
There's another bug related to this one IMHO:
[ 2202764 ] Taylor series of sqrt(1+x-y)

http://sourceforge.net/tracker/index.php?func=detail&aid=2202764&group_id=4933&atid=104933
[This is why i came here - i had a Taylor series bug that i think is
equivalent.]

HACK SOLUTION:
The bug is solved for me (maxima 5.10.0, debian-etch) by setting

  radexpand : false;

ANALYSIS:
However, i'm not convinced that this is a sufficient response to the
bug. i don't really know maxima well enough to know what would be
most reasonable.

(1) Set radexpand to false by default? And add more warnings in the
documentation (or do we expect users to learn mathematics elsewhere
than in software documentation?)

(2) Tell maxima not to rewrite sqrt(X)  as %i sqrt(-X) when domain:real ?

(3) Extend the "is X positive, negative or zero?" question to the
cases where it still needs to be asked but so far does not get asked,
i.e. when deciding whether

sqrt(X*Y) =  sqrt(-X)*sqrt(-Y)

or

sqrt(X*Y) = - sqrt(-X)*sqrt(-Y)  ?

(4) Decide not to invert the order of an expression inside sqrt( ... )
    if domain:real, except when it's sure that the expression is
positive?

----------------------------------------------------------------------

You can respond by visiting: 
https://sourceforge.net/tracker/?func=detail&atid=104933&aid=2184396&group_id=4933

[Maxima-bugs] [ maxima-Bugs-2184396 ] Wrong factorization of sqrt()

Computer Algebra System written in Common Lisp

[Maxima-bugs] [ maxima-Bugs-2184396 ] Wrong factorization of sqrt()