Where is the source for sage.matrix.matrix_modn_sparse.Matrix_modn_sparse.__eq__?

25 views
Skip to first unread message

Jeroen Demeyer

unread,
Oct 5, 2010, 9:29:11 PM10/5/10
to sage-devel
I would like to see the source code of
sage.matrix.matrix_modn_sparse.Matrix_modn_sparse.__eq__(), but I have
no idea where to look. The following doesn't work::

sage: L1 = matrix(GF(43), 3, 3, range(9), sparse=True)
sage: L1.__eq__??
Error getting source: arg is not a module, class, method, function,
traceback, frame, or code object
Type: method-wrapper
Base Class: <type 'method-wrapper'>
String Form: <method-wrapper '__eq__' of
sage.matrix.matrix_modn_sparse.Matrix_modn_sparse object at 0x5bb2dd0>
Namespace: Interactive
Docstring [source file open failed]:
x.__eq__(y) <==> x==y

François Bissey

unread,
Oct 6, 2010, 12:04:33 AM10/6/10
to sage-...@googlegroups.com
sage/matrix/matrix_modn_sparse.pyx ?

Jeroen Demeyer

unread,
Oct 6, 2010, 7:15:51 AM10/6/10
to sage-...@googlegroups.com
On 2010-10-06 02:04, Fran�ois Bissey wrote:
>> I would like to see the source code of
>> sage.matrix.matrix_modn_sparse.Matrix_modn_sparse.__eq__(), but I have
>> no idea where to look.

> sage/matrix/matrix_modn_sparse.pyx ?
Well, I obviously looked there, but didn't find it.

More precisely, I am looking for the location of the code which compares
a Matrix_modn_sparse with a Matrix_modn_dense. I.e. I want to see the
code for

sage: L = matrix(GF(43), 3, 3, range(9), sparse=True)
sage: R = matrix(GF(43), 3, 3, range(9))
sage: type(L)
<type 'sage.matrix.matrix_modn_sparse.Matrix_modn_sparse'>
sage: type(R)
<type 'sage.matrix.matrix_modn_dense.Matrix_modn_dense'>
sage: L==R # Where is this implemented?
True

I also tried the "trace()" command, but that doesn't seem to work well
with Cython code.

Jeroen.

François Bissey

unread,
Oct 6, 2010, 7:48:33 AM10/6/10
to sage-...@googlegroups.com
I was thinking that may be it was done by the __richcmp__, line 229-230:
def __richcmp__(matrix.Matrix self, right, int op): # always need for
mysterious reasons.
return self._richcmp(right, op)

Francois


Simon King

unread,
Oct 6, 2010, 8:12:12 AM10/6/10
to sage-devel
Hi Jeroen!

On Oct 6, 8:15 am, Jeroen Demeyer <jdeme...@cage.ugent.be> wrote:
> More precisely, I am looking for the location of the code which compares
> a Matrix_modn_sparse with a Matrix_modn_dense. I.e. I want to see the
> code for
>
> sage: L = matrix(GF(43), 3, 3, range(9), sparse=True)
> sage: R = matrix(GF(43), 3, 3, range(9))
> sage: type(L)
> <type 'sage.matrix.matrix_modn_sparse.Matrix_modn_sparse'>
> sage: type(R)
> <type 'sage.matrix.matrix_modn_dense.Matrix_modn_dense'>
> sage: L==R # Where is this implemented?
> True

If I am not mistaken, there is a generic __cmp__ (double underscore)
method at work, that will first send both arguments to a common
parent, by means of Sage's coercion model, and that then calls _cmp_
(single underscore).

So, if you ever want to implement comparison of elements of a parent
structure, I recall that you must *not* implement __cmp__ but _cmp_.

In this case, the common parent is the matrix space of *dense* 3x3
matrices over GF(43).

So, relevant is the _cmp_ method of the matrix R (since while
comparing L and R, L is sent to R's parent).

"R._cmp_?" shows that the comparison is forwarded to yet another
method (WHY??), namely _cmp_c_impl.

And I am afraid I can not find _cmp_c_impl.

So, it is indeed confusing.

Cheers,
Simon

John Cremona

unread,
Oct 6, 2010, 8:30:40 AM10/6/10
to sage-...@googlegroups.com

Indeed. Perhaps Jeroen's main point is that it is supposed to be
always possible to find out ow something is implemented using ??, but
in some cases like this that's not possible -- instead you do have to
know more about the structure of the Sage library.

A second point: if Simon is right then whenever a dense matrix is
compared to a sparse one, the sparse one is first converted to a dense
equivalent. That does not look very efficient to me, though it is
true that when both are of size mxn then all m*n entries of the dense
one might need to be looked at.

John

>
> Cheers,
> Simon
>
> --
> To post to this group, send an email to sage-...@googlegroups.com
> To unsubscribe from this group, send an email to sage-devel+...@googlegroups.com
> For more options, visit this group at http://groups.google.com/group/sage-devel
> URL: http://www.sagemath.org
>

Jeroen Demeyer

unread,
Oct 6, 2010, 8:51:20 AM10/6/10
to sage-...@googlegroups.com
On 2010-10-06 10:30, John Cremona wrote:
> Perhaps Jeroen's main point is that it is supposed to be
> always possible to find out ow something is implemented using ??, but
> in some cases like this that's not possible -- instead you do have to
> know more about the structure of the Sage library.
That wasn't really my main point (I actually want to find the source of
the == operator because I found a bug(*)). But of course you are right,
it would be nice if ?? would always give useful information (but that's
probably not easy to do).

Jeroen.


(*) I am working on improving (essentially rewriting) the interrupt
handling in Sage, see #9678. As an important part of this effort, I
have written code to test interrupt handling, see #10030. For the
moment, I have fixed some missing _sig_off's in #10061. I know there is
a missing _sig_off in the L == R code for matrices.

Jason Grout

unread,
Oct 6, 2010, 10:55:46 AM10/6/10
to sage-...@googlegroups.com
On 10/6/10 2:15 AM, Jeroen Demeyer wrote:
> On 2010-10-06 02:04, Fran�ois Bissey wrote:
>>> I would like to see the source code of
>>> sage.matrix.matrix_modn_sparse.Matrix_modn_sparse.__eq__(), but I have
>>> no idea where to look.
>
>> sage/matrix/matrix_modn_sparse.pyx ?
> Well, I obviously looked there, but didn't find it.
>
> More precisely, I am looking for the location of the code which compares
> a Matrix_modn_sparse with a Matrix_modn_dense. I.e. I want to see the
> code for
>
> sage: L = matrix(GF(43), 3, 3, range(9), sparse=True)
> sage: R = matrix(GF(43), 3, 3, range(9))
> sage: type(L)
> <type 'sage.matrix.matrix_modn_sparse.Matrix_modn_sparse'>
> sage: type(R)
> <type 'sage.matrix.matrix_modn_dense.Matrix_modn_dense'>
> sage: L==R # Where is this implemented?
> True
>


You can use the coercion model code to understand more about what
happens [1]:

sage: cm = sage.structure.element.get_coercion_model(); cm
<sage.structure.coerce.CoercionModel_cache_maps object at 0x101fb4770>
sage: sage: L = matrix(GF(43), 3, 3, range(9), sparse=True)
sage: sage: R = matrix(GF(43), 3, 3, range(9))
sage: cm.explain(L,R,operator.eq)
Equal but distinct parents.
Coercion on right operand via
Call morphism:
From: Full MatrixSpace of 3 by 3 dense matrices over Finite Field
of size 43
To: Full MatrixSpace of 3 by 3 sparse matrices over Finite
Field of size 43
Arithmetic performed after coercions.
Result lives in Full MatrixSpace of 3 by 3 sparse matrices over Finite
Field of size 43
Full MatrixSpace of 3 by 3 sparse matrices over Finite Field of size 43

This says that both are converted to sparse matrices before the
comparison happens. Looking at the code, we see that there is a
__richcmp__ (the cython equivalent of __cmp__ [2]) in the
matrix_modn_sparse.pyx file:

def __richcmp__(matrix.Matrix self, right, int op): # always need
for mysterious reasons.
return self._richcmp(right, op)

This is not defined in the matrix/ directory:

~/sage/devel/sage/sage/matrix% grep "def _richcmp[^_]" *.py*

~/sage/devel/sage/sage/matrix%

So we go up the chain and look in structure/ directory, and find it is
defined in devel/sage/sage/structure/element.pyx in the Element class,
which ends with:

return left._richcmp_c_impl(right, op)

and we also see a comment following the function:

####################################################################
# For a derived Cython class, you **must** put the following in
# your subclasses, in order for it to take advantage of the
# above generic comparison code. You must also define
# either _cmp_c_impl (if your subclass is totally ordered),
# _richcmp_c_impl (if your subclass is partially ordered), or both
# (if your class has both a total order and a partial order;
# then the total order will be available with cmp(), and the partial
# order will be available with the relation operators; in this case
# you must also define __cmp__ in your subclass).
# This is simply how Python works.


The _richcmp_c_impl is not defined in the matrix/ directory, so it must
also be inherited. Indeed, the _richcmp_c_impl method just below the
comment above is:

cdef _richcmp_c_impl(left, Element right, int op):
if (<Element>left)._richcmp_c_impl == Element._richcmp_c_impl and \
(<Element>left)._cmp_c_impl == Element._cmp_c_impl:
# Not implemented, try some basic defaults
if op == Py_EQ:
return left is right
elif op == Py_NE:
return left is not right
return left._rich_to_bool(op, left._cmp_c_impl(right))

so the _cmp_c_impl is called. So we need to go back down and see where
_cmp_c_impl. That's back in matrix/matrix_sparse.pyx:

cdef int _cmp_c_impl(self, Element right) except -2:
return cmp(self._dict(), right._dict())


So in the end, the dictionaries are compared, after the dense matrix is
coerced to a sparse matrix.


As for the "mysteriously needed" reason? Cython calls __richcmp__,
which needs to call the coercion framework _richcmp function. I guess
the mysterious part is that __richcmp__ is not inherited. Well,
according to the Python/C API docs [3], either all of __hash__, __cmp__,
and __richcmp__ are inherited, or none are. So I guess we need the
__richcmp__ function copied into the matrix/*.pyx files because __hash__
is redefined, so __richcmp__ isn't inherited.

At least, that's how I see it.

Thanks,

Jason

[1] http://www.sagemath.org/doc/reference/coercion.html, specifically
http://www.sagemath.org/doc/reference/coercion.html#basic-arithmetic-rules

[2]
http://docs.cython.org/src/userguide/special_methods.html#rich-comparisons

[3] see the thread
http://www.mail-archive.com/cytho...@codespeak.net/msg07374.html or
the docs
http://docs.python.org/c-api/typeobj.html?highlight=__hash__#tp_richcompare

Simon King

unread,
Oct 6, 2010, 11:21:17 AM10/6/10
to sage-devel
Hi John!

On Oct 6, 9:30 am, John Cremona <john.crem...@gmail.com> wrote:
> A second point: if Simon is right then whenever a dense matrix is
> compared to a sparse one, the sparse one is first converted to a dense
> equivalent.  That does not look very efficient to me,  though it is
> true that when both are of size mxn then all m*n entries of the dense
> one might need to be looked at.

I was slightly mistaken: There is no _cmp_ with single underscores.

In sage.structure.element, you will find the following comment.

####################################################################
# For a derived Cython class, you **must** put the following in
# your subclasses, in order for it to take advantage of the
# above generic comparison code. You must also define
# either _cmp_c_impl (if your subclass is totally ordered),
# _richcmp_c_impl (if your subclass is partially ordered), or both
# (if your class has both a total order and a partial order;
# then the total order will be available with cmp(), and the
partial
# order will be available with the relation operators; in this
case
# you must also define __cmp__ in your subclass).
# This is simply how Python works.
#
# For a *Python* class just define __cmp__ as always.
# But note that when this gets called you can assume that
# both inputs have identical parents.
#
# If your __cmp__ methods are not getting called, verify that the
# canonical_coercion(x,y) is not throwing errors.
#

####################################################################

So, the things are a little complicated, in the sense that the
behaviour for Cython and Python is different.

Moreover, if __cmp__ is implemented then it is called *after*
coercion.

Cheers,
Simon

David Joyner

unread,
Oct 6, 2010, 2:10:46 PM10/6/10
to sage-...@googlegroups.com
This is a really excellent response!

Do you think it should be added to the FAQ on the wiki ("I wanted to find
the code for a function using
sage: foo??
but couldn't. What is going on here?"),
or is it too technical?


On Wed, Oct 6, 2010 at 6:55 AM, Jason Grout <jason...@creativetrax.com> wrote:
> On 10/6/10 2:15 AM, Jeroen Demeyer wrote:
>>

John Cremona

unread,
Oct 6, 2010, 2:25:59 PM10/6/10
to sage-...@googlegroups.com
On Wed, Oct 6, 2010 at 3:10 PM, David Joyner <wdjo...@gmail.com> wrote:
> This is a really excellent response!
>
> Do you think it should be added to the FAQ on the wiki ("I wanted to find
> the code for a function using
> sage: foo??
> but couldn't. What is going on here?"),

+1

kcrisman

unread,
Oct 6, 2010, 2:46:03 PM10/6/10
to sage-devel


On Oct 6, 10:25 am, John Cremona <john.crem...@gmail.com> wrote:
> On Wed, Oct 6, 2010 at 3:10 PM, David Joyner <wdjoy...@gmail.com> wrote:
> > This is a really excellent response!
>
> > Do you think it should be added to the FAQ on the wiki ("I wanted to find
> > the code for a function using
> > sage: foo??
> > but couldn't. What is going on here?"),
>
> +1

Well, what should really be done is to make this sort of thing not
happen when we do ?? Jason, is that at all technically feasible, or
is it even more complicated than your coercion answer?

- kcrisman

Simon King

unread,
Oct 6, 2010, 7:10:51 PM10/6/10
to sage-devel
Hi all!
Yes, ?? should show the code. But, as far as I know, this is difficult
for special methods in Cython.

For example, if I am not mistaken, if X is an instance of a Cython
class, then X.__cmp__? (and similarly for other special methods) will
not even show you the docstring that was given to the method.

But even *if* ?? would show the source code: We must not forget that,
for example, X.__mul__ is a generic method inherited from a different
class, and there are good reasons for it (X being a ring element).

So, in these cases, showing the code would *not* show how the
functionality is ultimately implemented.

Or do you suggest that X.__mul__?? should also show the code of
X._mul_ ? I doubt that this would be feasible/reasonable.

I think
* ...? should show the documentation and ...?? show the code, even for
special methods in Cython, and also if a decorator is used; I don't
know if/how this can be achieved, though.
* If there is a generic method, like __mul__, then X.__mul__?? should
show the code of this generic method, but should not forward to the
code of a different method like X._mul_.
* The code of any generic method should provide extensive comments on
how this generic method should be implemented.

For example,
sage: RingElement.__mul__??
provides the comment
# Try fast pathway if they are both RingElements and the
parents match.
# (We know at least one of the arguments is a RingElement. So
if their
# types are *equal* (fast to check) then they are both
RingElements.
# Otherwise use the slower test via PY_TYPE_CHECK.)

That's a good thing (and better than Element.__mul__, that has no
documentation at all)!

But it could still be improved. I suggest that in addition it is
explicitly stated what a programmer should do, such as
# Multiplication of ring elements is supposed to be implemented in
a
# method called _mul_(self,other).
# That method returns the product of self and other. In that
method,
# self and other are guaranteed to belong to the same ring.
# Please do not overwrite the generic __mul__ method!

I know that this kind of information is available in some places in
the documentation. But some redundancy is helpful, IMO. Why not have
it stated in the code *and* in the docstring?

Cheers,
Simon
Reply all
Reply to author
Forward
Message has been deleted
0 new messages