Repository: sympy/sympy
Pull request: Fix Factors().as_expr(): don't multiply rational exponents
<!-- Your title above should be a short description of what
was changed. Do not include the issue number in the title. -->

#### References to other Issues or PRs
<!-- If this pull request fixes an issue, write "Fixes #NNNN" in that exact
format, e.g. "Fixes #1234". See
https://github.com/blog/1506-closing-issues-via-pull-requests . Please also
write a comment on that issue linking back to this pull request once it is
open. -->
Fixes #18059 

#### Brief description of what is fixed or changed
This commit fixes an invalid simplification in exprtools.py:

Before:
`Factors(sqrt(x**2)).as_expr() == x`
After:
`Factors(sqrt(x**2)).as_expr() == sqrt(x**2)`

This invalid simplification also led to a wrong result
by fu() (and hence trigsimp and simplify). Issue #18059.

Test cases are added.

#### Other comments


#### Release Notes

<!-- Write the release notes for this release below. See
https://github.com/sympy/sympy/wiki/Writing-Release-Notes for more information
on how to write release notes. The bot will check your release notes
automatically to see if they are formatted correctly. -->

<!-- BEGIN RELEASE NOTES -->
* core
  * Fixed a bug in `Factors().as_expr()` that led to simplifications that are not generally valid for complex numbers. (Exponents were always multiplied in power-of-a-power situations.)
<!-- END RELEASE NOTES -->


Problem statement:
Simplify of simple trig expression fails
trigsimp in various versions, including 1.5, incorrectly simplifies cos(x)+sqrt(sin(x)**2) as though it were cos(x)+sin(x) for general complex x. (Oddly it gets this right if x is real.)

Embarrassingly I found this by accident while writing sympy-based teaching material...

Hints from issue discussion:
I guess you mean this:
```julia
In [16]: cos(x) + sqrt(sin(x)**2)                                                                                                 
Out[16]: 
   _________         
  ╱    2             
╲╱  sin (x)  + cos(x)

In [17]: simplify(cos(x) + sqrt(sin(x)**2))                                                                                       
Out[17]: 
      ⎛    π⎞
√2⋅sin⎜x + ─⎟
      ⎝    4⎠
```
Which is incorrect if `sin(x)` is negative:
```julia
In [27]: (cos(x) + sqrt(sin(x)**2)).evalf(subs={x:-1})                                                                            
Out[27]: 1.38177329067604

In [28]: simplify(cos(x) + sqrt(sin(x)**2)).evalf(subs={x:-1})                                                                    
Out[28]: -0.301168678939757
```
For real x this works because the sqrt auto simplifies to abs before simplify is called:
```julia
In [18]: x = Symbol('x', real=True)                                                                                               

In [19]: simplify(cos(x) + sqrt(sin(x)**2))                                                                                       
Out[19]: cos(x) + │sin(x)│

In [20]: cos(x) + sqrt(sin(x)**2)                                                                                                 
Out[20]: cos(x) + │sin(x)│
```
Yes, that's the issue I mean.
`fu` and `trigsimp` return the same erroneous simplification. All three simplification functions end up in Fu's `TR10i()` and this is what it returns:
```
In [5]: from sympy.simplify.fu import *

In [6]: e = cos(x) + sqrt(sin(x)**2)

In [7]: TR10i(sqrt(sin(x)**2))
Out[7]: 
   _________
  ╱    2    
╲╱  sin (x) 

In [8]: TR10i(e)
Out[8]: 
      ⎛    π⎞
√2⋅sin⎜x + ─⎟
      ⎝    4⎠
```
The other `TR*` functions keep the `sqrt` around, it's only `TR10i` that mishandles it. (Or it's called with an expression outside its scope of application...)
I tracked down where the invalid simplification of `sqrt(x**2)` takes place or at least I think so:
`TR10i` calls `trig_split` (also in fu.py) where the line
https://github.com/sympy/sympy/blob/0d99c52566820e9a5bb72eaec575fce7c0df4782/sympy/simplify/fu.py#L1901
in essence applies `._as_expr()` to `Factors({sin(x)**2: S.Half})` which then returns `sin(x)`.

If I understand `Factors` (sympy.core.exprtools) correctly, its intent is to have an efficient internal representation of products and `.as_expr()` is supposed to reconstruct a standard expression from such a representation. But here's what it does to a general complex variable `x`:
```
In [21]: Factors(sqrt(x**2))
Out[21]: Factors({x**2: 1/2})
In [22]: _.as_expr()
Out[22]: x
```
It seems line 455 below
https://github.com/sympy/sympy/blob/0d99c52566820e9a5bb72eaec575fce7c0df4782/sympy/core/exprtools.py#L449-L458
unconditionally multiplies exponents if a power of a power is encountered. However this is not generally valid for non-integer exponents...

And line 457 does the same for other non-integer exponents:
```
In [23]: Factors((x**y)**z)
Out[23]: Factors({x**y: z})

In [24]: _.as_expr()
Out[24]:
 y⋅z
x
```

Commit message: Fix Factors().as_expr(): don't multiply rational exponents

This commit fixes an invalid simplification in exprtools.py:

Before:
`Factors(sqrt(x**2)).as_expr() == x`
After:
`Factors(sqrt(x**2)).as_expr() == sqrt(x**2)`

This invalid simplification also led to a wrong result
by fu() (and hence trigsimp and simplify). Issue #18059.

Test cases are added.

Patch to review:
diff --git a/sympy/core/exprtools.py b/sympy/core/exprtools.py
index 57fad7040760..a9db0d5c9839 100644
--- a/sympy/core/exprtools.py
+++ b/sympy/core/exprtools.py
@@ -448,14 +448,12 @@ def as_expr(self):  # Factors
         args = []
         for factor, exp in self.factors.items():
             if exp != 1:
-                b, e = factor.as_base_exp()
                 if isinstance(exp, int):
+                    b, e = factor.as_base_exp()
                     e = _keep_coeff(Integer(exp), e)
-                elif isinstance(exp, Rational):
-                    e = _keep_coeff(exp, e)
+                    args.append(b**e)
                 else:
-                    e *= exp
-                args.append(b**e)
+                    args.append(factor**exp)
             else:
                 args.append(factor)
         return Mul(*args)
diff --git a/sympy/core/tests/test_exprtools.py b/sympy/core/tests/test_exprtools.py
index d85ad06cf246..e482aad07e8f 100644
--- a/sympy/core/tests/test_exprtools.py
+++ b/sympy/core/tests/test_exprtools.py
@@ -27,6 +27,8 @@ def test_Factors():
     assert Factors({x: 2, y: 3, sin(x): 4}).as_expr() == x**2*y**3*sin(x)**4
     assert Factors(S.Infinity) == Factors({oo: 1})
     assert Factors(S.NegativeInfinity) == Factors({oo: 1, -1: 1})
+    # issue #18059:
+    assert Factors((x**2)**S.Half).as_expr() == (x**2)**S.Half
 
     a = Factors({x: 5, y: 3, z: 7})
     b = Factors({      y: 4, z: 3, t: 10})
diff --git a/sympy/simplify/tests/test_fu.py b/sympy/simplify/tests/test_fu.py
index 7d4a9c811726..dffeb1d289a9 100644
--- a/sympy/simplify/tests/test_fu.py
+++ b/sympy/simplify/tests/test_fu.py
@@ -276,6 +276,9 @@ def test_fu():
     expr = Mul(*[cos(2**i) for i in range(10)])
     assert fu(expr) == sin(1024)/(1024*sin(1))
 
+    # issue #18059:
+    assert fu(cos(x) + sqrt(sin(x)**2)) == cos(x) + sqrt(sin(x)**2)
+
 
 def test_objective():
     assert fu(sin(x)/cos(x), measure=lambda x: x.count_ops()) == \


Code review comments:
diff_hunk: @@ -448,14 +448,12 @@ def as_expr(self):  # Factors
         args = []
         for factor, exp in self.factors.items():
             if exp != 1:
-                b, e = factor.as_base_exp()
                 if isinstance(exp, int):
path: sympy/core/exprtools.py
text: @user1:
Could `exp` be a SymPy Integer?

@author:
Yes and in fact that seems to be the norm (at least for the few simple manual tests that I ran). In that case the general case applies (`args.append(factor**exp)`) and still returns the right result. I must admit though that I'm wondering a bit about that special-casing of the `int` case. Is that only to have it sympified? I'm also unsure about the `_keep_coeff`.
Any further insight is much appreciated.

@user2:
`_keep_coeff` is just there to keep 2-arg muls from evaluating (see docstring)

@user2:
From a quick look at Factors it seems that all values should be Basic type so the check for `int` should, perhaps, just be `Integer`.

@author:
What's strange though is that IIRC codecov marked the code for the `isinstance(exp, int)` case as being covered.
So both cases (int and Integer) should ideally end up using `_keep_coeff`?
I'll change it to `isinstance(exp, (int, Integer))` then.

@user2:
But the only place I can see an int potentially creeping in is at lines 360-361 where a default of 0 might be added to something that came from powers dict. If the raw 0 were changed to S.Zero then that would fix that (assuming that `as_powers_dict` might be returning an int. If the int case were fixed then problems like here will be avoided elsewhere.