Logged In: YES
user_id=895922
Originator: NO
As a workaround, you might try setting radexpand to false:
(%i1) assume(0 < x, x <= 1)$
(%i2) sqrt(a/(x-1)^2);
(%o2) sqrt(a)/(x-1)
(%i3) sqrt(a/(x-1)^2), radexpand : false;
(%o3) sqrt(a/(x-1)^2)
Dear Developers of Maxima,
I found that sqrt() returns the incorrect expression that has the sign
opposite to the ture expression when some certain argument is given to
sqrt(). Namely, sqrt() interprets incorrectly the database that is prepared
by assume().
Here is an demonstration program:
-------------------------------------------------------------------------------
/*
* wrong_sign_of_sqrt.maxima
*
* S.Adachi 2008/06/01
*/
display2d:false;
assume(x >= 0, x <= 1);
correct_result_1:sqrt((x-1)^2);
correct_result_2:sqrt(1/(x-1)^2);
correct_result_3:sqrt(a*(x-1)^2);
incorrect_result_1:sqrt(a/(x-1)^2);
incorrect_result_2:sqrt(a^2/(x-1)^2);
incorrect_result_3:sqrt(x^2/(x-1)^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_sign_of_sqrt.maxima)
batching #p/Volumes/HFS+2/home/adachi/work/299/wrong_sign_of_sqrt.maxima
(%i2) display2d : false
(%o2) false
(%i3) assume(x >= 0,x <= 1)
(%o3) [x >= 0,x <= 1]
(%i4) correct_result_1:sqrt((x-1)^2)
(%o4) 1-x
(%i5) correct_result_2:sqrt(1/(x-1)^2)
(%o5) 1/(1-x)
(%i6) correct_result_3:sqrt(a*(x-1)^2)
(%o6) sqrt(a)*(1-x)
(%i7) incorrect_result_1:sqrt(a/(x-1)^2)
(%o7) sqrt(a)/(x-1)
(%i8) incorrect_result_2:sqrt(a^2/(x-1)^2)
(%o8) abs(a)/(x-1)
(%i9) incorrect_result_3:sqrt(x^2/(x-1)^2)
(%o9) x/(x-1)
(%o10) "wrong_sign_of_sqrt.maxima"
-------------------------------------------------------------------------------
I wonder why sqrt() returns the wrong expression if sqrt((x-1)^2) appears
in the denominator of some fraction in the argument that is more complex
than some threshold (e.g. the numerator is not a simple number or something
like that).
I think that this is a very severe bug of sqrt()
and the database that is prepared by assume().
This bug puts many user programs to the state
producing many wrong results.
Please fix this severe bug.
Sincerely yours,
Satoshi Adachi
Logged In: YES
user_id=895922
Originator: NO
As a workaround, you might try setting radexpand to false:
(%i1) assume(0 < x, x <= 1)$
(%i2) sqrt(a/(x-1)^2);
(%o2) sqrt(a)/(x-1)
(%i3) sqrt(a/(x-1)^2), radexpand : false;
(%o3) sqrt(a/(x-1)^2)
Logged In: YES
user_id=1953419
Originator: YES
Thank you for your suggestion.
But, your suggestion just forbids Maxima to simplify certain expressions
in order not to produce wrong results.
Is someone trying to fix this bug now?
If not yet, I will read lisp source code in some future (maybe, several months
later) to try to fix the problem...
Logged In: YES
user_id=501686
Originator: NO
Assign category = lisp core / simplification.
Logged In: YES
user_id=28849
Originator: NO
FWIW, the code that handles this is simpexpt in simp.lisp. For the case of sqrt((x-1)^2), the e1 tag is done, and the second cond case is executed. The result is correct in this case.
However, for sqrt(1/(x-1)^2), the e1 tag is also done, but the first cond case is executed because the (noneg (cadr gr)) returns non-nil. In the former case the noneg call returns NIL, which is correct because x-1 is not negative.
I don't know why the two different calls to noneg return different values for x-1. NIL should be returned in each case because x <= 1.
There is one difference, though. The z arg to simpexpt is NIL for the case that works and T for the case that is incorrect.
Tracing dcompare shows that the results of comparisons with x
change in the middle of the computation.
(%i4) sqrt(a/(x-1)^2);
0: (DCOMPARE $X 1)
0: DCOMPARE returned $NZ
0: (DCOMPARE $X 1.0)
0: DCOMPARE returned $PNZ
0: (DCOMPARE $X 0)
0: DCOMPARE returned $PZ
0: (DCOMPARE $X 1)
0: DCOMPARE returned $POS
Similar behavior observed with bug 1419046 "sign bug".
Fixed in compar.lisp rev 1.46.
(%i2) assume(x >= 0, x <= 1);
(%o2) [x >= 0,x <= 1]
(%i3) correct_result_1:sqrt((x-1)^2);
(%o3) 1-x
(%i4) correct_result_2:sqrt(1/(x-1)^2);
(%o4) 1/(1-x)
(%i5) correct_result_3:sqrt(a*(x-1)^2);
(%o5) sqrt(a)*(1-x)
(%i6) incorrect_result_1:sqrt(a/(x-1)^2);
(%o6) sqrt(a)/(1-x)
(%i7) incorrect_result_2:sqrt(a^2/(x-1)^2);
(%o7) abs(a)/(1-x)
(%i8) incorrect_result_3:sqrt(x^2/(x-1)^2);
(%o8) x/(1-x)
Log in to post a comment.