On 26 February 2016 at 17:15, Aaron Meurer <asmeu...@gmail.com> wrote:
>> I traced this to here:
>>
>> In [1]: expr = x*sin(y)**2 + x*cos(y)**2 + 0.1*sin(y)**2 + 0.1*cos(y)**2
>>
>> In [2]: p = Poly(expr, x, cos(y), sin(y), domain='RR')
>>
>> In [3]: p
>> Out[3]: Poly(1.0*x*cos(y)**2 + 1.0*x*sin(y)**2 + 0.1*cos(y)**2 +
>> 0.1*sin(y)**2, x, cos(y), sin(y), domain='RR')
>>
>> In [4]: p.rep
>> Out[4]: DMP([[[1.0], [], [1.0, 0.0, 0.0]], [[0.1], [], [0.1, 0.0,
>> 0.0]]], RR, None)
>>
>> In [5]: p.rep.factor_list()
>> Out[5]:
>> (1.0,
>> [(DMP([[[0.1], [], [0.1, 0.0, 0.0]]], RR, None), 1),
>> (DMP([[[1.0]], [[0.1]]], RR, None), 1)])
>>
>> I think the bug is in there somehow but I have no idea what it means.
>> What is p.rep and how do you interpret it?
>
> It's the internal representation of the polynomial. It's a list of
> lists of lists of the coefficients (one level deep for each variable
> of the polynomial). factor_list() has a bug with floating point
> coefficients, it seems. Can you open an issue, if you haven't already?
I will do soon as I've almost fixed it. Actually I've found the bug in
dmp_factor_list in sympy/polys/factortools.py.
A simple example is
In [1]: x**2 - 0.01*y**2
Out[1]:
2 2
x - 0.01⋅y
In [2]: factor(x**2 - 0.01*y**2)
Out[2]: 1.0⋅(0.1⋅x - 0.01⋅y)⋅(0.1⋅x + 0.01⋅y)
In [3]: expand(factor(x**2 - 0.01*y**2))
Out[3]:
2 2
0.01⋅x - 0.0001⋅y
The bug is due to this:
https://github.com/sympy/sympy/blob/master/sympy/polys/factortools.py#L1272
The first pulls out a factor of 1/100 from the polynomial so that we get
(1/100) and (100*x**2 - y**2)
the second part is factored to
(10*x + y)*(10*x - y)
The factor of 100 is brought back in here:
https://github.com/sympy/sympy/blob/master/sympy/polys/factortools.py#L1304
However it is divided from BOTH factors to give
(10*x + y)/100 * (10*x - y)/100
So the result is 100x smaller than it should be. This will happen any
time you have non-integer valued floating point coefficients and there
are at least 2 polynomial factors in the expression.
I have an easy fix in dmp_factor_list:
$ git diff master
diff --git a/sympy/polys/factortools.py b/sympy/polys/factortools.py
index ebbabb8..ccecda6 100644
--- a/sympy/polys/factortools.py
+++ b/sympy/polys/factortools.py
@@ -1302,7 +1302,8 @@ def dmp_factor_list(f, u, K0):
coeff = coeff/denom
else:
for i, (f, k) in enumerate(factors):
- f = dmp_quo_ground(f, denom, u, K0)
+ if i == 0:
+ f = dmp_quo_ground(f, denom, u, K0)
f = dmp_convert(f, u, K0, K0_inexact)
factors[i] = (f, k)
I think I can see the same bug in dup_factor_list:
https://github.com/sympy/sympy/blob/master/sympy/polys/factortools.py#L1197
https://github.com/sympy/sympy/blob/master/sympy/polys/factortools.py#L1225
But I don't know how to trigger a code path that would go through
dup_factor_list. Under what circumstances would I hit that path?
--
Oscar
--
You received this message because you are subscribed to the Google Groups
"sympy" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to sympy+unsubscr...@googlegroups.com.
To post to this group, send email to sympy@googlegroups.com.
Visit this group at https://groups.google.com/group/sympy.
To view this discussion on the web visit
https://groups.google.com/d/msgid/sympy/CAHVvXxRufriYEVFd69udFzMwJPMaAepxQ8OhCUcmayksQEOdeA%40mail.gmail.com.
For more options, visit https://groups.google.com/d/optout.
