Cythonize matrix group elements

[tscrim]

They are used in a lot of places and should be cythonized. With my branch, I get ~10% speedup of my operations. Moreover, the class hierarchy for matrix group elements can and needs to be changed in order to do this. Specifically, I removed the MatrixGroupElement_base (the AffineGroupElement overrode all but 1 method) and GroupElementMixinLibGAP (only used by MatrixGroupElement_gap).

[tscrim]

While this new structure has some near-duplication by not using the calls to matrix(), it does away with a diamond problem and has fewer function calls overall by accessing the internal structure when necessary. Moreover, the original doctests for MatrixGroupElement_generic were not actually testing that class.

---

New commits:

​989b6ea	Adding changes from #19821.
​6358d53	Cythonized matrix_gp/group_element.py and simplified the class structure.

[chapoton]

many failing doctests, see patchbot

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​65ce465	Fixing modform_hecketriangle due to changes.
​1110926	Fixing last doctests.

[tscrim]

Fixed.

[vdelecroix]

For MatrixGroupElement_generic the __invert__ is very bad since it actually does a copy. If you have a matrix over ZZ its inverse belongs to matrices over QQ.

sage: m = matrix(3, [6, -1, -22, -5, 1, 18, -2, -4, 17])
sage: m.inverse().parent()
Full MatrixSpace of 3 by 3 dense matrices over Rational Field

Hence it would be better to use change_ring before the wrapping.

On the other hand there exists _invert_unit(). Compare

sage: m = matrix(3, [6, -1, -22, -5, 1, 18, -2, -4, 17])
sage: %timeit m.inverse().change_ring(ZZ)
10000 loops, best of 3: 53.7 µs per loop
sage: %timeit m._invert_unit()
100000 loops, best of 3: 6.56 µs per loop

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​3efe9b6	Merge branch 'develop' into public/groups/cythonize_matrix_group_element-19870
​29f7253	Special case for dense matrices over ZZ and making sure the inverse is in the base ring.

[tscrim]

At some point we do have to make a copy of the underlying matrix, and we can't modify self._matrix. It turns out that _invert_unit is only applicable to dense matrices over \ZZ, but it deserves a special case since checking if the base ring (in python) is less than 100 ns on my machine. I also made sure that the result lives in the base ring it should (which covers the sparse \ZZ case).

[vdelecroix]

Why not

try:
    M = M._invert_unit()
except AttributeError:
    M = ~M
    if M.base_ring() is not parent.base_ring():
        M = M.change_ring(parent.base_ring())

That would let the possibility for other rings/sparse matrices to implement an _invert_unit if appropriate. Though, if somebody changes the name of _invert_unit it will not be discovered. An alternative would be to implement a naive _invert_unit directly as a matrix method (basically following what you did).

[tscrim]

The exception handling is slower when it has to occur:

def foo(M):
    try:                                    
        return M._invert_unit()
    except AttributeError:
        return ~M

def bar(M):
    if M.base_ring() is ZZ and M.is_dense():
        return M._invert_unit()
    else:                 
        return ~M

sage: M = random_matrix(ZZ, 3, 3, sparse=True)
sage: M.det()
-1
sage: %timeit foo(M)
The slowest run took 15.06 times longer than the fastest. This could mean that an intermediate result is being cached 
1000 loops, best of 3: 643 µs per loop
sage: %timeit bar(M)
1000 loops, best of 3: 639 µs per loop
sage: M = random_matrix(QQ, 3, 3)
sage: M.det()
1
sage: %timeit foo(M)
The slowest run took 20.62 times longer than the fastest. This could mean that an intermediate result is being cached 
100000 loops, best of 3: 10.4 µs per loop
sage: %timeit bar(M)
The slowest run took 17.06 times longer than the fastest. This could mean that an intermediate result is being cached 
100000 loops, best of 3: 7.5 µs per loop

However, it is faster when no exception is generated:

sage: M = random_matrix(ZZ, 3, 3)
sage: M.det()
1
sage: %timeit foo(M)
The slowest run took 16.19 times longer than the fastest. This could mean that an intermediate result is being cached 
100000 loops, best of 3: 4.39 µs per loop
sage: %timeit bar(M)
The slowest run took 11.94 times longer than the fastest. This could mean that an intermediate result is being cached 
100000 loops, best of 3: 4.61 µs per loop

Yet, I don't think that is enough of a speedup to justify the other cases (in particular, the \QQ example).

[tscrim]

Actually, perhaps directly implementing a default _inverse_unit is best as it might be useful for matrices over other non-fields. Thoughts?

[vdelecroix]

I like the idea of a _inverse_unit but I would actually be in favor of _inverse_unit_. Moreover, such method makes sense for any ring, not only matrices. We just run into the same kind of trouble while implementing any (multiplicative) subgroup of units. The trivial example being

sage: parent(~ZZ(1))
Rational Field

[tscrim]

I agree that _inverse_unit_ is better since it is (more like) a special function. The quick solution would be to put a default for this in the Rings category or in Element which just does return self.parent()(~self). Or do you think we should do something a little more comprehensive on a separate ticket?

[jdemeyer]

1. I think you should add a .pxd file such that other Cython modules can cimport this if needed.

2. Use the standard ​copyright template.

3. Remove __cmp__ = _cmp_ twice (it doesn't work and it's not needed).

4. You don't actually use this: from sage.structure.element import have_same_parent.

5. There is not much point in making matrix() a cpdef method: you can just directly access the _matrix attribute from Cython.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​5333b0f

Implementing reviewer comments.

[tscrim]

Replying to jdemeyer:

1. I think you should add a .pxd file such that other Cython modules can cimport this if needed.

2. Use the standard ​copyright template.

3. Remove __cmp__ = _cmp_ twice (it doesn't work and it's not needed).

4. You don't actually use this: from sage.structure.element import have_same_parent.

5. There is not much point in making matrix() a cpdef method: you can just directly access the _matrix attribute from Cython.

All done. I also made _matrix a public attribute since doing access that way was fastest way to access the defining matrix by my testing.

[chapoton]

ok, good enough for me. And patchbots are green.

In a later ticker, we should overload the is_one method.

[vbraun]

Merge conflict, try the next beta...

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​df18b13

Merge branch 'public/groups/cythonize_matrix_group_element-19870' of trac.sagemath.org:sage into public/groups/cythonize_matrix_group_element-19870

[tscrim]

Trivial rebase because git forgot that the file changed names (even though I had done it with git mv).

[jdemeyer]

In the .pxd, you should not repeat methods coming from base classes. So the _cmp_ and _mul_ should not be there.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​ddd203e

Fixing group_element.pxd file by not redefining inherited methods.

[tscrim]

Replying to jdemeyer:

In the .pxd, you should not repeat methods coming from base classes. So the _cmp_ and _mul_ should not be there.

Done.

[tscrim]

Thank you.