RR, coercion of numpy.float32() and clang - coercion experts needed

85 views
Skip to first unread message

Dima Pasechnik

unread,
Apr 20, 2017, 5:10:29 PM4/20/17
to sage-devel
polynomials with numpy coefficients strike again, this time if you build Sage with clang:

sage: import numpy as np
sage: x=polygen(RR)        # this is problematic
sage: y=polygen(RDF)      # this works fine
sage: np.float('1.5')*x
1.50000000000000*x
sage: np.float('1.5')*y
1.5*x
sage: np.float32('1.5')*x
/usr/home/dima/Sage/sage/src/bin/sage-ipython:1: RuntimeWarning: invalid value encountered in multiply
  #!/usr/bin/env python
1.50000000000000*x
sage: np.float64('1.5')*x
/usr/home/dima/Sage/sage/src/bin/sage-ipython:1: RuntimeWarning: invalid value encountered in multiply
  #!/usr/bin/env python
1.50000000000000*x
sage: np.float64('1.5')*y
1.5*x
sage: np.float32('1.5')*y
1.5*x


Finding out where exactly in numpy this warning happens is not easy (one can set up a numpy interrupt handler, but because it happens in compiled code, error message is almost meaningless)
but it looks as if it might be a coercion problem.
Any ideas where to look?

Marc Mezzarobba

unread,
Apr 21, 2017, 4:11:37 AM4/21/17
to sage-...@googlegroups.com
Dima Pasechnik wrote:
> but it looks as if it might be a coercion problem.
> Any ideas where to look?

Not really, but it does look like the common parent
discovered by the coercion system is incorrect:

sage: import numpy as np
sage: a = np.float('1.5')
sage: b = np.float32('1.5')

sage: get_coercion_model().common_parent(b, polygen(RR))
Univariate Polynomial Ring in x over Real Field with 53 bits of
precision

sage: RR.coerce(a)
1.50000000000000

sage: RR.coerce(b)
...
TypeError: no canonical coercion from <type 'numpy.float32'> to Real
Field with 53 bits of precision

sage: get_coercion_model().common_parent(b, RR)
<type 'numpy.float32'>

--
Marc

Dima Pasechnik

unread,
Apr 21, 2017, 5:10:12 AM4/21/17
to sage-devel, ma...@mezzarobba.net
This does not look like this is the root cause of the problem.
Namely, note that  also
 
sage: a128 = np.float128('1.5') 
sage: RR.coerce(a128)
# throws the same as above TypeError, but
sage: a128*x
# does not print the numpy warning.

what is also somewhat puzzling is 
that while 
sage: np.float32('1.5')*x
# does print the numpy warning

sage: x*np.float32('1.5')
# just works (i.e. does not print the numpy warning)

So this is a subtle combination of an apparent bug in numpy (or in clang) with 
some Sage coercion weirdness.
Apparently numpy's floats are coerced into RDF for the purpose of
dealing with polynomials, see py_scalar_to_element in structure/coerce.pyx, but apparently this only happens
for 
 x*np.float32('1.5'),  and not for np.float32('1.5')*x ?

Dima



--
Marc

David Roe

unread,
Apr 21, 2017, 10:46:03 AM4/21/17
to sage-devel, ma...@mezzarobba.net
Quick guess (don't have time to look at this soon): the issue could be in the action discovery code, which may try creating elements of numpy.float32 using weird inputs. This could cause a difference between right and left multiplication. One way to check if this is the source is to try converting each of the polynomials from RR['x'].some_elements() into numpy.float32 and see if you get the same error. 

--
You received this message because you are subscribed to the Google Groups "sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+unsubscribe@googlegroups.com.
To post to this group, send email to sage-...@googlegroups.com.
Visit this group at https://groups.google.com/group/sage-devel.
For more options, visit https://groups.google.com/d/optout.

Dima Pasechnik

unread,
Apr 21, 2017, 11:29:26 AM4/21/17
to sage-devel, ma...@mezzarobba.net


On Friday, April 21, 2017 at 3:46:03 PM UTC+1, David Roe wrote:
Quick guess (don't have time to look at this soon): the issue could be in the action discovery code, which may try creating elements of numpy.float32 using weird inputs. This could cause a difference between right and left multiplication. One way to check if this is the source is to try converting each of the polynomials from RR['x'].some_elements() into numpy.float32 and see if you get the same error. 

no, this does not give anything that looks like the error in question (no errors on degree 0 terms,
"ValueError: setting an array element with a sequence." on degree >0 terms)

 

To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.

David Roe

unread,
Apr 21, 2017, 9:10:48 PM4/21/17
to sage-devel, ma...@mezzarobba.net
On Fri, Apr 21, 2017 at 11:29 AM, Dima Pasechnik <dim...@gmail.com> wrote:


On Friday, April 21, 2017 at 3:46:03 PM UTC+1, David Roe wrote:
Quick guess (don't have time to look at this soon): the issue could be in the action discovery code, which may try creating elements of numpy.float32 using weird inputs. This could cause a difference between right and left multiplication. One way to check if this is the source is to try converting each of the polynomials from RR['x'].some_elements() into numpy.float32 and see if you get the same error. 

no, this does not give anything that looks like the error in question (no errors on degree 0 terms,
"ValueError: setting an array element with a sequence." on degree >0 terms)

Sorry, I don't have any other ideas, and can't play around with it easily since the standard Sage build isn't producing these errors.
David
To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+unsubscribe@googlegroups.com.

Dima Pasechnik

unread,
Apr 22, 2017, 4:22:21 AM4/22/17
to sage-devel
it is not hard to build with clang on OSX, see #12426.

and this issue is causing the majority of doctest failures in this setting, as well as on clang/FreeBSD.(#22679).
(where the building is not as smooth yet).

I suppose we will have to merge these clang-enabling changes one way or another soon.

Dima Pasechnik

unread,
Apr 22, 2017, 5:01:52 AM4/22/17
to sage-devel, ma...@mezzarobba.net, roed...@gmail.com


On Saturday, April 22, 2017 at 2:10:48 AM UTC+1, David Roe wrote:


On Fri, Apr 21, 2017 at 11:29 AM, Dima Pasechnik <dim...@gmail.com> wrote:


On Friday, April 21, 2017 at 3:46:03 PM UTC+1, David Roe wrote:
Quick guess (don't have time to look at this soon): the issue could be in the action discovery code, which may try creating elements of numpy.float32 using weird inputs. This could cause a difference between right and left multiplication. One way to check if this is the source is to try converting each of the polynomials from RR['x'].some_elements() into numpy.float32 and see if you get the same error. 

no, this does not give anything that looks like the error in question (no errors on degree 0 terms,
"ValueError: setting an array element with a sequence." on degree >0 terms)

Sorry, I don't have any other ideas, and can't play around with it easily since the standard Sage build isn't producing these errors.
David


Is there a "normal" way to trace all the coersion/conversion steps for a given Sage expression?
If not, there should be one, IMHO (otherwise it's always a huge guessing game for this sort of issues)

Volker Braun

unread,
Apr 22, 2017, 5:47:48 AM4/22/17
to sage-devel, ma...@mezzarobba.net, roed...@gmail.com
There is this:

sage: import numpy as np
sage: cm = sage.structure.element.coercion_model
sage: cm.explain(polygen(RR), np.float32('1.5'), operator.mul)
Action discovered.
    Right scalar multiplication by Real Double Field on Univariate Polynomial Ring in x over Real Field with 53 bits of precision
    with precomposition on right by Native morphism:
      From: Set of Python objects of type 'numpy.float32'
      To:   Real Double Field
Result lives in Univariate Polynomial Ring in x over Real Field with 53 bits of precision
Univariate Polynomial Ring in x over Real Field with 53 bits of precision


Daniel Krenn

unread,
Apr 22, 2017, 5:50:31 AM4/22/17
to sage-...@googlegroups.com
On 2017-04-22 11:01, Dima Pasechnik wrote:
> Is there a "normal" way to trace all the coersion/conversion steps for a
> given Sage expression?
> If not, there should be one, IMHO (otherwise it's always a huge guessing
> game for this sort of issues)

What helps is .record_exceptions and .coercion_traceback:

http://doc.sagemath.org/html/en/reference/structure/sage/structure/element.html#sage.structure.element.coercion_traceback

Daniel

Dima Pasechnik

unread,
Apr 22, 2017, 6:22:59 PM4/22/17
to sage-devel, ma...@mezzarobba.net, roed...@gmail.com


On Saturday, April 22, 2017 at 10:47:48 AM UTC+1, Volker Braun wrote:
There is this:

sage: import numpy as np
sage: cm = sage.structure.element.coercion_model
sage: cm.explain(polygen(RR), np.float32('1.5'), operator.mul)
Action discovered.
    Right scalar multiplication by Real Double Field on Univariate Polynomial Ring in x over Real Field with 53 bits of precision
    with precomposition on right by Native morphism:
      From: Set of Python objects of type 'numpy.float32'
      To:   Real Double Field
Result lives in Univariate Polynomial Ring in x over Real Field with 53 bits of precision
Univariate Polynomial Ring in x over Real Field with 53 bits of precision

Thanks! It seems that cm.explain() lies in this case, i.e., the action is actually computed in a different way!!!
Indeed:

sage: import numpy as np
sage: x=polygen(RealField())
sage: np.float32('1.5')*x
/usr/home/dima/Sage/sage/src/bin/sage-ipython:1: RuntimeWarning: invalid value encountered in multiply
  #!/usr/bin/env python
1.50000000000000*x
sage: cm = sage.structure.element.get_coercion_model()
sage: cm.explain(np.float32('1.5'), x, operator.mul)
Action discovered.
    Left scalar multiplication by Real Double Field on Univariate Polynomial Ring in x over Real Field with 53 bits of precision
    with precomposition on left by Native morphism:
      From: Set of Python objects of type 'numpy.float32'
      To:   Real Double Field
Result lives in Univariate Polynomial Ring in x over Real Field with 53 bits of precision
Univariate Polynomial Ring in x over Real Field with 53 bits of precision
sage: RDF(np.float32('1.5'))*x
1.50000000000000*x

that is, cm.explain tells me that   np.float32('1.5')*x is computed by doing
RDF(np.float32('1.5'))*x

But this is not the case, as RDF(np.float32('1.5'))*x does not produce any numpy warnings, whereas
np.float32('1.5')*x does.

Then I went digging for  the place which is supposed to do apply RDF(), and tried to break  it:

--- a/src/sage/structure/coerce.pyx
+++ b/src/sage/structure/coerce.pyx
@@ -231,7 +231,8 @@ cpdef py_scalar_to_element(x):
             return Integer(x)
         elif isinstance(x, numpy.floating):
             from sage.rings.real_double import RDF
-            return RDF(x)
+           # return RDF(x)
+            return RDF(42)
         elif isinstance(x, numpy.complexfloating):
             from sage.rings.complex_double import CDF
             return CDF(x)

but to no avail, 1.5 does not get replaced by 42.0 in the answer.

It seems that this part of coersion is broken in the sense that it does something rogue,
not prescribed by the framework---and it has nothing to do with clang (I also tried this "patch" in
a usual configuration).
Reply all
Reply to author
Forward
0 new messages