#22924 closed enhancement (fixed)

cleaning of linbox for dense integer matrices

Reported by: vdelecroix Owned by:
Priority: major Milestone: sage-8.0
Component: packages: standard Keywords: linbox
Cc: cpernet, dlucas Merged in:
Authors: Vincent Delecroix Reviewers: Clément Pernet
Report Upstream: N/A Work issues:
Branch: 0b082c2 (Commits) Commit: 0b082c2d7aac83d91253a450cd81e7accf13d76b
Dependencies: Stopgaps:

Description (last modified by vdelecroix)

Various cleaning to the linbox interface (sage.libs.linbox) and Sage dense integer matrices (sage.matrix.matrix_integer_dense).

  • implement a more direct linbox/flint interface in sage/libs/linbox/linbox_flint_interface.pyx (supersedes part of linbox-sage.C inside linbox)
  • remove all mpz pointers from Matrix_integer_dense (that is the attributes bint _initialized_mpz, mpz_t * _entrie and mpz_t ** _rows) as well as the associated constructors/destructors (int _init_mpz(self), int _init_mpz_impl(self), int _init_linbox(self), void _dealloc_mpz(self))
  • fix the following error
    sage: matrix(ZZ, 2).det(algorithm='padic')
    Traceback (most recent call last):
    ...
    ImportError: No module named matrix_integer_dense_hnf
    
  • add more examples and tests to the documentation

(see #22872 for the linbox task ticket)

Change History (49)

comment:1 Changed 22 months ago by vdelecroix

  • Branch set to u/vdelecroix/22924
  • Commit set to 38d9da303a72512d2cab94ef0d8219fd03de12aa

New commits:

cca8edffix headers
38d9da3Cleaning linbox interface of Matrix_integer_dense

comment:2 Changed 22 months ago by git

  • Commit changed from 38d9da303a72512d2cab94ef0d8219fd03de12aa to 86d8af3dd66ed8c1847452d555953eee7a00f91c

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

86d8af3Cleaning linbox interface of Matrix_integer_dense

comment:3 Changed 22 months ago by vdelecroix

  • Description modified (diff)
  • Status changed from new to needs_review

comment:4 Changed 22 months ago by vdelecroix

  • Description modified (diff)

comment:5 Changed 22 months ago by vdelecroix

  • Component changed from PLEASE CHANGE to packages: standard

comment:6 Changed 22 months ago by git

  • Commit changed from 86d8af3dd66ed8c1847452d555953eee7a00f91c to 4cae968e969b661488ef07621d29303251732726

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

cccbe30Cleaning linbox interface of Matrix_integer_dense
4cae96822924: remove (instead of deprecate) code from linbox.pyx

comment:7 Changed 22 months ago by git

  • Commit changed from 4cae968e969b661488ef07621d29303251732726 to 09f62806bc83b5049d90ea481de8abffe61712f0

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

09f628022924: get rid of spies

comment:8 Changed 22 months ago by git

  • Commit changed from 09f62806bc83b5049d90ea481de8abffe61712f0 to 69a140e1d299cc6971ad75e868dd9418083bbc0b

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

48cd83bfix headers
3f5ec83Cleaning linbox interface of Matrix_integer_dense
ef07ed522924: remove (instead of deprecate) code from linbox.pyx
69a140e22924: get rid of spies

comment:9 Changed 22 months ago by vdelecroix

(rebased on 8.0.b4)

comment:10 follow-up: Changed 22 months ago by cpernet

Great. Let's start with a few remarks:

  1. Why do need the dependency on #22886 ? That ticket is stalled and is therefore also blocking this one.
  2. While we are changing a lot of stuff, I wouldn't mind getting rid of these ugly poly_linbox method, with a string parameter typ that determines whether to compute minpoly or charpoly. We should rather only have the 2 methods _minpoly_linbox_ and _charpoly_linbox_. There is not much duplicated code between these two, especially now after this clean-up.

comment:11 Changed 22 months ago by cpernet

  • Reviewers set to Clément Pernet

comment:12 in reply to: ↑ 10 ; follow-up: Changed 22 months ago by vdelecroix

Replying to cpernet:

Great.

Thanks for having a look!

Let's start with a few remarks:

  1. Why do need the dependency on #22886 ? That ticket is stalled and is therefore also blocking this one.

I think that Cython compilation fails otherwise (I am trying again right now on top of 8.0.beta4... will take some time since I need to recompile the whole Sage after playing with #22493).

  1. While we are changing a lot of stuff, I wouldn't mind getting rid of these ugly poly_linbox method, with a string parameter typ that determines whether to compute minpoly or charpoly. We should rather only have the 2 methods _minpoly_linbox_ and _charpoly_linbox_. There is not much duplicated code between these two, especially now after this clean-up.

The ticket only intends to clean matrix_integer_dense.pyx. The two functions you mention are actually removed from matrix/matrix_integer_dense.pyx but kept in libs/linbox/linbox_flint_interface.pyx (but note that the argument there is an integer). I can get rid of the following function in the interface

# set p to the minpoly or the charpoly of A
cdef inline void linbox_fmpz_mat_poly(fmpz_poly_t p, fmpz_mat_t A, int minpoly)

Is that what you want?

comment:13 in reply to: ↑ 12 Changed 22 months ago by cpernet

Replying to vdelecroix:

I can get rid of the following function in the interface

# set p to the minpoly or the charpoly of A
cdef inline void linbox_fmpz_mat_poly(fmpz_poly_t p, fmpz_mat_t A, int minpoly)

Is that what you want?

Ok sorry for the confusion. My remark indeed applies to the new libs/linbox_flint_interface.pyx. Sure, it should definitely not be exposed through the new interface as it is an internal trick to factor out code. Then I also suggest to avoid this trick and simply write the two methods for the sake of clarity.

comment:14 Changed 22 months ago by git

  • Commit changed from 69a140e1d299cc6971ad75e868dd9418083bbc0b to 637fd872a78a7aec0781998a281d4dbbbc0fd7af

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

55bdc0aCleaning linbox interface of Matrix_integer_dense
c8b95de22924: remove (instead of deprecate) code from linbox.pyx
fcafa6022924: get rid of spies
637fd8722924: cleaning in linbox_flint_interface.pyx

comment:15 Changed 22 months ago by vdelecroix

  • Dependencies #22886 deleted

No dependency on #22886 anymore (it is needed for cylinbox but it is fine for this ticket).

I cleaned the flint interface a bit more than what you asked for.

comment:16 Changed 22 months ago by vdelecroix

I was wrong (and you were right). There is still def _poly_linbox(self, typ, var='x'): as a method of Matrix_integer_dense.

comment:17 follow-up: Changed 22 months ago by vdelecroix

By the way, I am adding tests and LinBox has trouble with the characteristic polynomial of the 0 x 0 matrix

sage: matrix([]).charpoly('y', 'linbox')
Traceback (most recent call last):
...
MemoryError: failed to allocate 2305843008945258512 bytes

Does it qualify as a bug? For now, I will add a special case in the interface LinBox/flint.

EDIT: and it hangs forever for minpoly.

Last edited 22 months ago by vdelecroix (previous) (diff)

comment:18 follow-up: Changed 22 months ago by vdelecroix

For minpoly it is written

        .. NOTE::

           Linbox charpoly disabled on 64-bit machines, since it hangs
           in many cases.

However, I do not see where this actually happens in the code. Do you know what this is?

comment:19 Changed 22 months ago by dlucas

  • Cc dlucas added

comment:20 Changed 22 months ago by vdelecroix

Un petit timing

sage: m = random_matrix(ZZ, 50)
sage: %timeit m._clear_cache(); p = m.charpoly(algorithm='linbox')
100 loops, best of 3: 16.8 ms per loop
sage: %timeit m._clear_cache(); p = m.charpoly(algorithm='flint')
10 loops, best of 3: 44.2 ms per loop
sage: %timeit m._clear_cache(); p = m.charpoly(algorithm='generic')
10 loops, best of 3: 152 ms per loop
Last edited 22 months ago by vdelecroix (previous) (diff)

comment:21 Changed 22 months ago by git

  • Commit changed from 637fd872a78a7aec0781998a281d4dbbbc0fd7af to fde1353a5c22bd3ff3fe50c9df982264b30d505e

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

fde135322924: clean minpoly/charpoly in Matrix_integer_dense

comment:22 Changed 22 months ago by vdelecroix

All right: Matrix_integer_dense cleaned and interface fixed on corner cases with commit fde1353.

comment:23 Changed 22 months ago by vdelecroix

And patchbot is happy (excepted some timeout I am not responsible for).

comment:24 in reply to: ↑ 18 Changed 22 months ago by cpernet

Replying to vdelecroix:

For minpoly it is written

        .. NOTE::

           Linbox charpoly disabled on 64-bit machines, since it hangs
           in many cases.

However, I do not see where this actually happens in the code. Do you know what this is?

It's a super old comment. Some archeology on trac tickets should tell where it comes from but it's safe to remove it.

comment:25 follow-up: Changed 22 months ago by cpernet

There are a couple of warnings that one could avoid:

[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:313:49: local variable 'res' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:389:46: local variable 'res' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:394:42: local variable 'res' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:398:42: local variable 'res' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:498:34: local variable 'g' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:532:32: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:534:16: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:562:30: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:564:16: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:717:33: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:736:33: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:756:35: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:1110:39: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:1127:39: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:1150:39: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:1167:39: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:1184:39: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:1217:35: local variable 'r' referenced before assignment
[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:1218:59: local variable 'r' referenced before assignment
[

Basically, the variables are passed as "out" arguments, and cython complains that they are not initialized. It is safe to force their value to 0 at declaration (e.g. int res = 0;) and this would silence out the warning.

comment:26 in reply to: ↑ 17 Changed 22 months ago by cpernet

Replying to vdelecroix:

By the way, I am adding tests and LinBox has trouble with the characteristic polynomial of the 0 x 0 matrix

sage: matrix([]).charpoly('y', 'linbox')
Traceback (most recent call last):
...
MemoryError: failed to allocate 2305843008945258512 bytes

Does it qualify as a bug? For now, I will add a special case in the interface LinBox/flint.

EDIT: and it hangs forever for minpoly.

This is definitely a bug. Reported upstream: https://github.com/linbox-team/linbox/issues/51

comment:27 follow-up: Changed 22 months ago by cpernet

In linbox_flint_interface.pyx lines 84 and 157, you seem to use the method size of givvector to get the degree of the corresponding polynomial. While this is very likely to be correct in the current code, the specification do not require that the size of the vector storing the polynomial is ajusted to its degree. One should rather use Givaro::Poly1Dom::degree() instead.

comment:28 in reply to: ↑ 25 ; follow-ups: Changed 22 months ago by vdelecroix

Replying to cpernet:

There are a couple of warnings that one could avoid:

[sagelib-8.0.beta4] warning: sage/rings/finite_rings/element_givaro.pyx:313:49: local variable 'res' referenced before assignment
...

Basically, the variables are passed as "out" arguments, and cython complains that they are not initialized. It is safe to force their value to 0 at declaration (e.g. int res = 0;) and this would silence out the warning.

Yes, but these warnings have nothing to do with this ticket. Feel free to open another one to deal with rings/finite_rings (and put me in cc if you do).

comment:29 in reply to: ↑ 27 Changed 22 months ago by vdelecroix

Replying to cpernet:

In linbox_flint_interface.pyx lines 84 and 157, you seem to use the method size of givvector to get the degree of the corresponding polynomial. While this is very likely to be correct in the current code, the specification do not require that the size of the vector storing the polynomial is ajusted to its degree. One should rather use Givaro::Poly1Dom::degree() instead.

In LinBox code there are two occurrences of v.size() in examples/charpoly.C for a const Polynomial &v (lines 47 and 57). Idem in examples/minpoly.C (line 36). Are these also wrong?

How am I supposed to use degree? I was able to do it 7 lines in Cython... it looks extremely complicated to get such elementary quantity.

cdef GivaroIntegerRing ZZ
cdef GivaroDegree deg
cdef LinBoxIntegerPolynomialRing * R
R = new LinBoxIntegerPolynomialRing(ZZ)
R.degree(deg, q)
del R
cdef size_t length = deg.value() + 1

comment:30 Changed 22 months ago by git

  • Commit changed from fde1353a5c22bd3ff3fe50c9df982264b30d505e to 6b45893cfec10b3bc20cc44a64233f66e0cf0d8a

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

a75ec9122924: better FIXME message
6b4589322924: compute length from degree

comment:31 follow-up: Changed 22 months ago by vdelecroix

I implemented your remarks. Though, I do not think that 6b45893 is reasonable.

comment:32 in reply to: ↑ 31 Changed 22 months ago by cpernet

Replying to vdelecroix:

I implemented your remarks. Though, I do not think that 6b45893 is reasonable.

Ok this seems indeed quite tedious. After browsing through the current use of size() for degree in LinBox?, I would say it's rather safe to assume that the vector has been resized and therefore that the size method would return the degree. I was probably being too pedantic, sorry. I suggest to drop 6b45893.

Last edited 22 months ago by cpernet (previous) (diff)

comment:33 in reply to: ↑ 28 Changed 22 months ago by cpernet

Replying to vdelecroix:

Yes, but these warnings have nothing to do with this ticket. Feel free to open another one to deal with rings/finite_rings (and put me in cc if you do).

My bad, I thought they were introduced in your new code. Will open a new ticket

comment:34 Changed 22 months ago by git

  • Commit changed from 6b45893cfec10b3bc20cc44a64233f66e0cf0d8a to a75ec910414f17772072ea13d57fb11117877c16

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

comment:35 follow-up: Changed 22 months ago by cpernet

Also, I see some new code not directly related to the subject of the ticket (e.g. in _rational_kernel_flint, rational_kernel_iml), but this is rather ok for me as it is documentation improvement on routines related to the ticket).

Last edited 22 months ago by cpernet (previous) (diff)

comment:36 follow-up: Changed 22 months ago by vdelecroix

Replying to cpernet:

Replying to vdelecroix:

I implemented your remarks. Though, I do not think that 6b45893 is reasonable.

Ok this seems indeed quite tedious. After browsing through the current use of size() for degree in LinBox?, I would say it's rather safe to assume that the vector has been resized and therefore that the size method would return the degree. I was probably being too pedantic, sorry. I suggest to drop 6b45893.

I dropped it.

It would be fine for me to use degree however it is horribly complicated. Why is there a class Degree that is a trivial wrapper around a int64_t? Why could we not use degree directly as a method of a polynomial?

comment:37 in reply to: ↑ 35 Changed 22 months ago by vdelecroix

Replying to cpernet:

Also, I see some new code not directly related to the subject of the ticket (e.g. in _rational_kernel_flint, rational_kernel_iml), but this is rather ok for me as it is documentation improvement on routines related to the ticket).

It is. _rational_kernel_iml now uses fmpz_mat_to_mpz_array that was introduced in the ticket.

comment:38 in reply to: ↑ 28 Changed 22 months ago by cpernet

Replying to vdelecroix:

Feel free to open another one to deal with rings/finite_rings (and put me in cc if you do).

Done in #22966.

Last edited 22 months ago by cpernet (previous) (diff)

comment:39 in reply to: ↑ 36 ; follow-up: Changed 22 months ago by cpernet

Let me try an explanation for this mess (which I agree should be simplified):

Replying to vdelecroix:

Why is there a class Degree that is a trivial wrapper around a int64_t?

I would have said: so as to deal with the degree of the 0 polynomial which is -1 and not size()-1=0

Why could we not use degree directly as a method of a polynomial?

I assume because in Givaro (and the whole LinBox ecosystem), the container (here givvector does not know about the semantic). It is the domain of computation (here GivPoly1Dom who does).

comment:40 in reply to: ↑ 39 ; follow-up: Changed 22 months ago by vdelecroix

Replying to cpernet:

Let me try an explanation for this mess (which I agree should be simplified):

Replying to vdelecroix:

Why is there a class Degree that is a trivial wrapper around a int64_t?

I would have said: so as to deal with the degree of the 0 polynomial which is -1 and not size()-1=0

In LinBox you have size(0) = 1!? Why would you store the zero at all?

Why could we not use degree directly as a method of a polynomial?

I assume because in Givaro (and the whole LinBox ecosystem), the container (here givvector does not know about the semantic). It is the domain of computation (here GivPoly1Dom who does).

I see. This way is indeed more flexible.

comment:41 Changed 22 months ago by vdelecroix

Tell me if there are other improvements to be done on this ticket.

comment:42 in reply to: ↑ 40 Changed 22 months ago by cpernet

Replying to vdelecroix:

I would have said: so as to deal with the degree of the 0 polynomial which is -1 and not size()-1=0

In LinBox you have size(0) = 1!? Why would you store the zero at all?

My mistake. Givaro actually resizes the 0 polynomial to an empty vector, so size-1 is always equal to the degree.

After a discussion with J-G. Dumas on the topic here are a few clarifications:

  • the degree can be accessed in two ways:
    1. either through `Poly1Dom<..>::degree which silently resize the vector when there are leading zeros.
    2. direcrly by P.size()-1 when the programmer knows what she does and does not want to pay the price of a resize.

So it's fine to do it with size()-1 as in this ticket, since we know the polynomials handed by LinBox? have correct size in this case.

comment:43 Changed 22 months ago by cpernet

  • Status changed from needs_review to positive_review

Ok, the ticket is now fine for me. -> positive review.

Last edited 22 months ago by cpernet (previous) (diff)

comment:44 Changed 22 months ago by cpernet

  • Status changed from positive_review to needs_work

Sorry, two minor documentation typos to be fixed in src/sage/libs/linbox/linbox_flint_interface.pyx:

  • line 10 "linbox-Sage.C" -> "linbox-sage.C"
  • line 128 "preformed" -> "performed"

Back in needs work.

comment:45 Changed 22 months ago by git

  • Commit changed from a75ec910414f17772072ea13d57fb11117877c16 to 0b082c2d7aac83d91253a450cd81e7accf13d76b

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

0b082c222924: typos

comment:46 Changed 22 months ago by vdelecroix

  • Status changed from needs_work to needs_review

comment:47 Changed 22 months ago by cpernet

  • Status changed from needs_review to positive_review

Fine then, back to positive review.

comment:48 Changed 22 months ago by vdelecroix

Thanks for the review!

comment:49 Changed 21 months ago by vbraun

  • Branch changed from u/vdelecroix/22924 to 0b082c2d7aac83d91253a450cd81e7accf13d76b
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.