use LinBox as native matrix representation for dense matrices over GF(p)

malb commented 15 years ago

Copying to and from LinBox uses up precious RAM and the point of fast linear algebra is to deal with large matrices. We should consider switching to LinBox as the native representation of matrices over GF(p)

Without Patch

sage: A = random_matrix(GF(97),2000,2000)
sage: %time A*A
CPU times: user 9.66 s, sys: 0.12 s, total: 9.77 s
Wall time: 9.82 s

With Patch

sage: A = random_matrix(GF(97),2000,2000)
sage: %time A*A
CPU times: user 1.32 s, sys: 0.00 s, total: 1.32 s
Wall time: 1.35 s

Magma

> A:=RandomMatrix(GF(97),2000,2000);
> time C:=A*A;                      
Time: 1.560

Install http://sage.math.washington.edu/home/malb/spkgs/linbox-1.1.6.p5.spkg
Apply attachment: trac_4260-linbox_default.patch
Apply attachment: trac_4260-dense_ctypes_template.patch
Apply attachment: trac_4260-matrix-modn-docs.patch
Apply attachment: trac_4260_more_doctests.patch
Apply attachment: trac_4260_echelonformdomain.patch
Apply attachment: trac_4260-minor_fixes.patch
Apply attachment: trac_4260_bugfix.patch

CC: @simon-king-jena @rbeezer @sagetrac-drkirkby

Component: linear algebra

Keywords: linbox, sd32, sd34

Author: Burcin Erocal, Martin Albrecht, Rob Beezer

Reviewer: Burcin Erocal, Simon King, Martin Albrecht, Jeroen Demeyer

Merged: sage-4.8.alpha3

Issue created by migration from https://trac.sagemath.org/ticket/4260

ClementPernet commented 15 years ago

comment:1

I will work on it as a coding sprint at SD10.

burcin commented 12 years ago

Author: Burcin Erocal

burcin commented 12 years ago

comment:3

I finally rebased the patch from SD16. The template class in the patch contains the updates made to the modndense class since then, like changes to the `sig*functions. Apparently the modn_dense class representation now allows permuting the rows by changing pointers in the_matrixarray. We can't allow that if we want to pass the_entries` to linbox, so I skipped those changes.

Sage builds with the attached patches, and you can construct matrices. However, there are lots of bugs, some linbox wrappers are still stubs, etc. Expect crashes and wrong results.

With the patch applied, I get a crash with the following:

sage: a = matrix(GF(97),3,4,range(12))
sage: a.echelonize()
*** glibc detected *** python: free(): invalid next size (fast): 0x000000000270b370 ***
======= Backtrace: =========
<snip>

AFAICT, the new cython code is an exact copy of the wrapper function in linbox-sage.C. Here is what valgrind says:

==3026== Invalid write of size 8
==3026==    at 0x39E49EF1: T.4552 (ffpack_ludivine.inl:420)
==3026==    by 0x39E49AA0: T.4552 (ffpack_ludivine.inl:486)
==3026==    by 0x39E4ABBF: __pyx_pf_4sage_6matrix_24matrix_modn_dense_double_26M
atrix_modn_dense_template_20_echelonize_linbox(_object*, _object*) (ffpack.h:113
2)
==3026==    by 0x4E74082: PyObject_Call (abstract.c:2492)
==3026==    by 0x39E2CA8A: __pyx_pf_4sage_6matrix_24matrix_modn_dense_double_26M
atrix_modn_dense_template_19echelonize(_object*, _object*, _object*) (matrix_mod
n_dense_double.cpp:8738)
==3026==    by 0x4F17FF9: PyEval_EvalFrameEx (ceval.c:3706)
==3026==    by 0x4F19CDC: PyEval_EvalCodeEx (ceval.c:2968)
==3026==    by 0x4F19DB1: PyEval_EvalCode (ceval.c:522)
==3026==    by 0x4F19083: PyEval_EvalFrameEx (ceval.c:4401)
==3026==    by 0x4F19CDC: PyEval_EvalCodeEx (ceval.c:2968)
<snip lots more Py* lines>
==3026==  Address 0x6ca16e8 is 0 bytes after a block of size 24 alloc'd
==3026==    at 0x4C267CE: malloc (vg_replace_malloc.c:236)
==3026==    by 0x39E4AA5A: __pyx_pf_4sage_6matrix_24matrix_modn_dense_double_26Matrix_modn_dense_template_20_echelonize_linbox(_object*, _object*) (memory.h:32)
==3026==    by 0x4E74082: PyObject_Call (abstract.c:2492)
==3026==    by 0x39E2CA8A: __pyx_pf_4sage_6matrix_24matrix_modn_dense_double_26Matrix_modn_dense_template_19echelonize(_object*, _object*, _object*) (matrix_modn_dense_double.cpp:8738)
==3026==    by 0x4F17FF9: PyEval_EvalFrameEx (ceval.c:3706)
==3026==    by 0x4F19CDC: PyEval_EvalCodeEx (ceval.c:2968)
==3026==    by 0x4F19DB1: PyEval_EvalCode (ceval.c:522)
==3026==    by 0x4F19083: PyEval_EvalFrameEx (ceval.c:4401)
==3026==    by 0x4F19CDC: PyEval_EvalCodeEx (ceval.c:2968)
==3026==    by 0x4F18074: PyEval_EvalFrameEx (ceval.c:3802)
<snip lots of Py* lines>

I'd appreciate any pointers about the problem above, though I don't know if I'll have the time to come back to this before the bug days in August (when I presume Martin will take over?).

1d7ec08f-60ae-4512-91a6-8324c06eab9f commented 12 years ago

comment:4

These are the files in sage/matrix with failures:

        sage -t  devel/sage-main/sage/matrix/matrix_cyclo_dense.pyx # 22 doctests failed
        sage -t  devel/sage-main/sage/matrix/strassen.pyx # 2 doctests failed
        sage -t  devel/sage-main/sage/matrix/matrix0.pyx # 2 doctests failed
        sage -t  devel/sage-main/sage/matrix/matrix_integer_dense.pyx # 5 doctests failed
        sage -t  devel/sage-main/sage/matrix/matrix_space.py # 1 doctests failed
        sage -t  devel/sage-main/sage/matrix/matrix_window_modn_dense.pyx # 1 doctests failed
        sage -t  devel/sage-main/sage/matrix/matrix_modn_sparse.pyx # 1 doctests failed
        sage -t  devel/sage-main/sage/matrix/matrix_integer_dense_saturation.py # 0 doctests failed
        sage -t  devel/sage-main/sage/matrix/matrix_rational_dense.pyx # 44 doctests failed
        sage -t  devel/sage-main/sage/matrix/matrix2.pyx # Time out
        sage -t  devel/sage-main/sage/matrix/matrix_modn_dense.pyx # Time out
        sage -t  devel/sage-main/sage/matrix/matrix_modn_dense_template.pxi # Time out

malb commented 12 years ago

Attachment: trac_4260-linbox_default.patch.gz

make matrix space constructor use the new classes

malb commented 12 years ago

Description changed:

--- 
+++ 
@@ -1 +1,25 @@
 Copying to and from LinBox uses up precious RAM and the point of fast linear algebra is to deal with large matrices. We should consider switching to LinBox as the native representation of matrices over GF(p)
+
+**Without Patch**
+
+```
+!python
+sage: A = random_matrix(GF(97),2000,2000)
+sage: %time A*A
+CPU times: user 9.66 s, sys: 0.12 s, total: 9.77 s
+Wall time: 9.82 s
+```
+
+**With Patch**
+
+```
+!python
+sage: A = random_matrix(GF(97),2000,2000)
+sage: %time A*A
+CPU times: user 1.32 s, sys: 0.00 s, total: 1.32 s
+Wall time: 1.35 s
+```
+
+* **Install** http://sage.math.washington.edu/home/malb/spkgs/linbox-1.1.6.p4.spkg
+* **Apply** [attachment: trac_4260-linbox_default.patch](https://github.com/sagemath/sage-prod/files/10642185/trac_4260-linbox_default.patch.gz)
+* **Apply** [attachment: trac_4260-dense_ctypes_template.patch](https://github.com/sagemath/sage-prod/files/10642186/trac_4260-dense_ctypes_template.patch.gz)

malb commented 12 years ago

comment:5

I fixed a few issues and segfaults but the thing is far from done. However, one can probably do higher level stuff now, i.e. it shouldn't crash that much any more.

We need a new LinBox SPKG because Modular<float> didn't have a NonZeroRandIter which is needed by the charpoly code. LinBox 1.1.7 fixes this issue but I tried unsuccessfully to upgrade to 1.1.7 for like 10 hours (cf. #11718).

malb commented 12 years ago

Description changed:

--- 
+++ 
@@ -2,8 +2,7 @@

 **Without Patch**

-```
-!python
+```python
 sage: A = random_matrix(GF(97),2000,2000)
 sage: %time A*A
 CPU times: user 9.66 s, sys: 0.12 s, total: 9.77 s
@@ -12,8 +11,7 @@

 **With Patch**

-```
-!python
+```python
 sage: A = random_matrix(GF(97),2000,2000)
 sage: %time A*A
 CPU times: user 1.32 s, sys: 0.00 s, total: 1.32 s

malb commented 12 years ago

comment:7

Doctest failures with most recent patch on sage.math:

        sage -t  -long -force_lib devel/sage/doc/de/tutorial/tour_advanced.rst # 2 doctests failed
        sage -t  -long -force_lib devel/sage/doc/en/tutorial/tour_advanced.rst # 2 doctests failed
        sage -t  -long -force_lib devel/sage/doc/en/bordeaux_2008/modular_forms_and_hecke_operators.rst # 1 doctests failed
        sage -t  -long -force_lib devel/sage/doc/en/bordeaux_2008/elliptic_curves.rst # 4 doctests failed
        sage -t  -long -force_lib devel/sage/doc/fr/tutorial/tour_advanced.rst # 2 doctests failed
        sage -t  -long -force_lib devel/sage/doc/ru/tutorial/tour_advanced.rst # 2 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modsym/heilbronn.pyx # 2 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modsym/tests.py # 1 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modsym/subspace.py # 9 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modsym/space.py # 18 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modform/eisenstein_submodule.py # 3 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modform/tests.py # 1 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modform/constructor.py # 3 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modform/space.py # 8 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modform/ambient.py # 4 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modform/cuspidal_submodule.py # 6 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modsym/ambient.py # 4 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/modform/element.py # 11 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/hecke/element.py # 1 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/hecke/hecke_operator.py # 1 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/hecke/module.py # 3 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/abvar/homology.py # 3 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/hecke/submodule.py # 3 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/abvar/torsion_subgroup.py # 4 doctests failed
        sage -t  -long -force_lib devel/sage/sage/modular/abvar/abvar.py # 4 doctests failed
        sage -t  -long -force_lib devel/sage/sage/matrix/matrix_cyclo_dense.pyx # 8 doctests failed
        sage -t  -long -force_lib devel/sage/sage/matrix/matrix2.pyx # 1 doctests failed
        sage -t  -long -force_lib devel/sage/sage/tests/cmdline.py # 1 doctests failed
        sage -t  -long -force_lib devel/sage/sage/combinat/symmetric_group_representations.py # 1 doctests failed
        sage -t  -long -force_lib devel/sage/sage/schemes/elliptic_curves/padics.py # 29 doctests failed
        sage -t  -long -force_lib devel/sage/sage/schemes/elliptic_curves/padic_lseries.py # 6 doctests failed
        sage -t  -long -force_lib devel/sage/sage/schemes/elliptic_curves/ell_modular_symbols.py # 2 doctests failed
        sage -t  -long -force_lib devel/sage/sage/schemes/generic/toric_chow_group.py # 16 doctests failed
        sage -t  -long -force_lib devel/sage/sage/schemes/elliptic_curves/ell_rational_field.py # 1 doctests failed
        sage -t  -long -force_lib devel/sage/sage/schemes/elliptic_curves/sha_tate.py # 10 doctests failed

malb commented 12 years ago

Changed author from Burcin Erocal to Burcin Erocal, Martin Albrecht

malb commented 12 years ago

Description changed:

--- 
+++ 
@@ -18,6 +18,14 @@
 Wall time: 1.35 s

+Magma + + +> A:=RandomMatrix(GF(97),2000,2000); +> time C:=A*A; +Time: 1.560 + +

Install http://sage.math.washington.edu/home/malb/spkgs/linbox-1.1.6.p4.spkg
Apply attachment: trac_4260-linbox_default.patch
Apply attachment: trac_4260-dense_ctypes_template.patch

malb commented 12 years ago

comment:8

I've fixed all the easy stuff which brings the doctest failures down to:

sage -t  -long devel/sage/doc/en/bordeaux_2008/elliptic_curves.rst # 4 doctests failed
sage -t  -long devel/sage/sage/modular/modsym/subspace.py # 9 doctests failed
sage -t  -long devel/sage/sage/modular/modsym/space.py # 18 doctests failed
sage -t  -long devel/sage/sage/modular/modform/eisenstein_submodule.py # 3 doctests failed
sage -t  -long devel/sage/sage/modular/modform/constructor.py # 3 doctests failed
sage -t  -long devel/sage/sage/modular/modform/space.py # 8 doctests failed
sage -t  -long devel/sage/sage/modular/modform/ambient.py # 4 doctests failed
sage -t  -long devel/sage/sage/modular/hecke/element.py # 1 doctests failed
sage -t  -long devel/sage/sage/modular/hecke/hecke_operator.py # 1 doctests failed
sage -t  -long devel/sage/sage/modular/hecke/module.py # 3 doctests failed
sage -t  -long devel/sage/sage/modular/abvar/homology.py # 3 doctests failed
sage -t  -long devel/sage/sage/modular/hecke/submodule.py # 3 doctests failed
sage -t  -long devel/sage/sage/modular/abvar/torsion_subgroup.py # 4 doctests failed
sage -t  -long devel/sage/sage/modular/abvar/abvar.py # 4 doctests failed
sage -t  -long devel/sage/sage/combinat/symmetric_group_representations.py # 1 doctests failed
sage -t  -long devel/sage/sage/schemes/elliptic_curves/padics.py # 29 doctests failed
sage -t  -long devel/sage/sage/schemes/elliptic_curves/padic_lseries.py # 6 doctests failed
sage -t  -long devel/sage/sage/schemes/elliptic_curves/ell_modular_symbols.py # 2 doctests failed
sage -t  -long devel/sage/sage/schemes/generic/toric_chow_group.py # 16 doctests failed
sage -t  -long devel/sage/sage/schemes/elliptic_curves/ell_rational_field.py # 1 doctests failed
sage -t  -long devel/sage/sage/schemes/elliptic_curves/sha_tate.py # 10 doctests failed

many of which seem to be caused by a small set of bugs.

malb commented 12 years ago

comment:9

Here's where we are at on sage.math:

sage -t  devel/sage/doc/en/bordeaux_2008/elliptic_curves.rst # 4 doctests failed
sage -t  devel/sage/sage/modular/modsym/subspace.py # 9 doctests failed
sage -t  devel/sage/sage/modular/modsym/space.py # 12 doctests failed
sage -t  devel/sage/sage/modular/modform/eisenstein_submodule.py # 1 doctests failed
sage -t  devel/sage/sage/modular/modform/space.py # 7 doctests failed
sage -t  devel/sage/sage/modular/modform/constructor.py # 1 doctests failed
sage -t  devel/sage/sage/modular/modform/ambient.py # 4 doctests failed
sage -t  devel/sage/sage/modular/hecke/element.py # 1 doctests failed
sage -t  devel/sage/sage/modular/hecke/hecke_operator.py # 1 doctests failed
sage -t  devel/sage/sage/modular/hecke/module.py # 3 doctests failed
sage -t  devel/sage/sage/modular/abvar/homology.py # 3 doctests failed
sage -t  devel/sage/sage/modular/abvar/torsion_subgroup.py # 4 doctests failed
sage -t  devel/sage/sage/modular/hecke/submodule.py # 3 doctests failed
sage -t  devel/sage/sage/modular/abvar/abvar.py # 4 doctests failed
sage -t  devel/sage/sage/structure/sage_object.pyx # 1 doctests failed
sage -t  devel/sage/sage/combinat/symmetric_group_representations.py # 1 doctests failed
sage -t  devel/sage/sage/schemes/elliptic_curves/padics.py # 29 doctests failed
sage -t  devel/sage/sage/schemes/elliptic_curves/padic_lseries.py # 6 doctests failed
sage -t  devel/sage/sage/schemes/elliptic_curves/ell_modular_symbols.py # 2 doctests failed
sage -t  devel/sage/sage/schemes/generic/toric_chow_group.py # 16 doctests failed
sage -t  devel/sage/sage/schemes/elliptic_curves/sha_tate.py # 10 doctests failed
sage -t  devel/sage/sage/schemes/elliptic_curves/ell_rational_field.py # 1 doctests failed

williamstein commented 12 years ago

Work Issues: sd32

malb commented 12 years ago

comment:11

With the updated patch we are down to:

sage -t  -long devel/sage/sage/modular/modsym/heilbronn.pyx # 2 doctests failed
sage -t  -long devel/sage/sage/modular/abvar/homology.py # 3 doctests failed

However, there also seems to be a doctest failure in matrix2.pyx which is not that easily reproduced.

malb commented 12 years ago

comment:12

Now all doctests should pass!

williamstein commented 12 years ago

Changed keywords from linbox to linbox, sd32

williamstein commented 12 years ago

Changed work issues from sd32 to none

malb commented 12 years ago

Work Issues: extend documentation

malb commented 12 years ago

Changed author from Burcin Erocal, Martin Albrecht to Burcin Erocal, Martin Albrecht, Rob Beezer

malb commented 12 years ago

Description changed:

--- 
+++ 
@@ -29,3 +29,4 @@
 * **Install** http://sage.math.washington.edu/home/malb/spkgs/linbox-1.1.6.p4.spkg
 * **Apply** [attachment: trac_4260-linbox_default.patch](https://github.com/sagemath/sage-prod/files/10642185/trac_4260-linbox_default.patch.gz)
 * **Apply** [attachment: trac_4260-dense_ctypes_template.patch](https://github.com/sagemath/sage-prod/files/10642186/trac_4260-dense_ctypes_template.patch.gz)
+* **Apply** [attachment: trac_4096-matrix-modn-docs.patch](https://github.com/sagemath/sage/files/ticket4260/trac_4096-matrix-modn-docs.patch.gz)

malb commented 12 years ago

Attachment: trac_4260-dense_ctypes_template.patch.gz

add templated classes for float and double representation

malb commented 12 years ago

comment:15

I adapted the crossover from float to double. Around 2¹¹ Modular<float>} is really slow because there are not enough bits left to let ATLAS do it's magic, i.e., too many modular reductions. On my computer using Modular<float> up to 2⁸ seems like a good choice. On sage.math this choice isn't too bad (but not optimal). Multiplying two 1,000 x 1,000 matrices over GF(p) (2nd column) which is smaller than 2ⁱ (1st column) and the time it takes:

 2          3 0.22000
 3          7 0.24000
 4         13 0.24000
 5         31 0.25000
 6         61 0.26000
 7        127 0.26000
 8        251 0.62000
 9        509 0.38000 <=== using Modular<double> now
10       1021 0.38000
11       2039 0.39000
12       4093 0.39000
13       8191 0.40000
14      16381 0.41000
15      32749 0.41000
16      65521 0.42000
17     131071 0.43000
18     262139 0.43000
19     524287 0.44000
20    1048573 0.44000
21    2097143 0.45000
22    4194301 0.66000
23    8388593 1.91000 <=== Generic matrices

simon-king-jena commented 12 years ago

comment:16

I found that the time for computing echelon form became worse:

sage-4.6.2

sage: MS = MatrixSpace(GF(101),2000,2000)
sage: %time A = MS.random_element()
CPU times: user 0.17 s, sys: 0.03 s, total: 0.20 s
Wall time: 0.20 s
sage: B = MS.random_element()
sage: %time C = A*B
CPU times: user 8.35 s, sys: 0.07 s, total: 8.42 s
Wall time: 8.45 s
sage: %time A.echelonize()
CPU times: user 1.22 s, sys: 0.06 s, total: 1.28 s
Wall time: 1.38 s

sage-4.7.2.alpha2 with the patches and spkg from here:

sage: MS = MatrixSpace(GF(101),2000,2000)
sage: %time A = MS.random_element()
CPU times: user 0.19 s, sys: 0.03 s, total: 0.22 s
Wall time: 0.22 s
sage: B = MS.random_element()
sage:  %time C = A*B
CPU times: user 1.16 s, sys: 0.02 s, total: 1.17 s
Wall time: 1.22 s
sage: %time A.echelonize()
CPU times: user 1.87 s, sys: 0.00 s, total: 1.87 s
Wall time: 1.92 s

malb commented 12 years ago

Attachment: trac_4260-matrix-modn-docs.patch.gz

Attachment: trac_4260_more_doctests.patch.gz

malb commented 12 years ago

Description changed:

--- 
+++ 
@@ -29,4 +29,5 @@
 * **Install** http://sage.math.washington.edu/home/malb/spkgs/linbox-1.1.6.p4.spkg
 * **Apply** [attachment: trac_4260-linbox_default.patch](https://github.com/sagemath/sage-prod/files/10642185/trac_4260-linbox_default.patch.gz)
 * **Apply** [attachment: trac_4260-dense_ctypes_template.patch](https://github.com/sagemath/sage-prod/files/10642186/trac_4260-dense_ctypes_template.patch.gz)
-* **Apply** [attachment: trac_4096-matrix-modn-docs.patch](https://github.com/sagemath/sage/files/ticket4260/trac_4096-matrix-modn-docs.patch.gz)
+* **Apply** [attachment: trac_4260-matrix-modn-docs.patch](https://github.com/sagemath/sage-prod/files/10642187/trac_4260-matrix-modn-docs.patch.gz)
+* **Apply** [attachment: trac_4260-more_doctests.patch](https://github.com/sagemath/sage/files/ticket4260/trac_4260-more_doctests.patch.gz)

malb commented 12 years ago

comment:17

renamed Rob's patch to fix ticket number
- (hopefully) added doctests for every single function
- Simon, can you try again after setting MAX_MODULUS in sage.matrix.matrix_modn_dense_float to 2⁶? This forces the use of doubles for GF(101) which might be more efficient. Also, how fast is A.echelonize('gauss') for you on that benchmark?

malb commented 12 years ago

Description changed:

--- 
+++ 
@@ -30,4 +30,4 @@
 * **Apply** [attachment: trac_4260-linbox_default.patch](https://github.com/sagemath/sage-prod/files/10642185/trac_4260-linbox_default.patch.gz)
 * **Apply** [attachment: trac_4260-dense_ctypes_template.patch](https://github.com/sagemath/sage-prod/files/10642186/trac_4260-dense_ctypes_template.patch.gz)
 * **Apply** [attachment: trac_4260-matrix-modn-docs.patch](https://github.com/sagemath/sage-prod/files/10642187/trac_4260-matrix-modn-docs.patch.gz)
-* **Apply** [attachment: trac_4260-more_doctests.patch](https://github.com/sagemath/sage/files/ticket4260/trac_4260-more_doctests.patch.gz)
+* **Apply** [attachment: trac_4260_more_doctests.patch](https://github.com/sagemath/sage-prod/files/10642188/trac_4260_more_doctests.patch.gz)

simon-king-jena commented 12 years ago

comment:19

Replying to @malb:

Simon, can you try again after setting MAX_MODULUS in sage.matrix.matrix_modn_dense_float to 2⁶? This forces the use of doubles for GF(101) which might be more efficient.

It isn't:

sage: sage.matrix.matrix_modn_dense_float.MAX_MODULUS = 2^6
sage: MS = MatrixSpace(GF(101),2000,2000)
sage: %time A = MS.random_element()
CPU times: user 0.21 s, sys: 0.01 s, total: 0.22 s
Wall time: 0.22 s
sage: B = MS.random_element()
sage: %time C = A*B
CPU times: user 1.88 s, sys: 0.04 s, total: 1.92 s
Wall time: 1.93 s
sage: %time A.echelonize()
CPU times: user 2.65 s, sys: 0.00 s, total: 2.65 s
Wall time: 2.69 s
sage: type(A)
<type 'sage.matrix.matrix_modn_dense_double.Matrix_modn_dense_double'>

Also, how fast is A.echelonize('gauss') for you on that benchmark?

You mean "how slow", I suppose:

sage: A = MS.random_element()
sage: %time A.echelonize('gauss')
CPU times: user 41.53 s, sys: 0.10 s, total: 41.63 s
Wall time: 41.75 s

malb commented 12 years ago

comment:20

Okay, so both the old code and this patch call LinBox but with the patch it's slower (I can reproduce your timing difference). Hence, we'll have to check what LinBox in the old version ends up doing vs. what it is doing now.

malb commented 12 years ago

Changed work issues from extend documentation to improve echelonize

ClementPernet commented 12 years ago

comment:21

A word about the regression (I'm copying my reply to malb on linbox-devel)

The new code (that I wrote )

size_t r = FFPACK::ReducedRowEchelonForm(F, nrows, ncols, matrix, ncols, P,Q);

calls the actual RowEchelon elimination in FFPACK, which transforms A into its redrowechelon form E and the transformation matrix U (both matrices being magically stored inplace in A)

It is slower than the older code sage-4.6 using linbox-1.1.6:

int rank = EF.rowReducedEchelon(E, A);

The latter computes the redrowechlon (actually the trans of the redcolechelon), but no transformation matrix. This saves roughly 50% of the total number of arithmetic ops (1n^3 rather than 2n^3), and explains the regression.

Switching back to the old way should fix the regression (for a quick fix). And I still need to add the feature of not computing the transform at the level of FFPACK, since I expect some timing improvements over the old version in linbox 1.1.6.

malb commented 12 years ago

comment:22

Hi Clement,

thanks for explaining. I always forget about the transformation matrix (e.g., that in fact M4RI is even faster than Magma than previously thought, because we always compute the transformation matrix and yet we are faster :)). I'll try to "switch back". Btw. is there a way to construct the right matrices without copying?

PS: I didn't get your reply on [linbox-devel] btw.

malb commented 12 years ago

Description changed:

--- 
+++ 
@@ -31,3 +31,4 @@
 * **Apply** [attachment: trac_4260-dense_ctypes_template.patch](https://github.com/sagemath/sage-prod/files/10642186/trac_4260-dense_ctypes_template.patch.gz)
 * **Apply** [attachment: trac_4260-matrix-modn-docs.patch](https://github.com/sagemath/sage-prod/files/10642187/trac_4260-matrix-modn-docs.patch.gz)
 * **Apply** [attachment: trac_4260_more_doctests.patch](https://github.com/sagemath/sage-prod/files/10642188/trac_4260_more_doctests.patch.gz)
+* **Apply** [attachment: trac_4260_echelonformdomain.patch](https://github.com/sagemath/sage-prod/files/10642190/trac_4260_echelonformdomain.patch.gz)

malb commented 12 years ago

comment:23

The additional patch makes us use the old EchelonFormDomain interface which is twice as fast (as Clément explained). Simon, can you review this ticket?

burcin commented 12 years ago

Attachment: trac_4260_echelonformdomain.patch.gz

Attachment: trac_4260-minor_fixes.patch.gz

minor fixes

burcin commented 12 years ago

comment:24

attachment: trac_4260-minor_fixes.patch

makes some cosmetic changes and
fixes a possible memory leak if allocation of these matrices fail.

I read through the patches and the resulting code. All looks good to me. Please switch this to positive review if my patch is ok.

Thanks everyone for finally finishing this off.

burcin commented 12 years ago

Changed keywords from linbox, sd32 to linbox, sd32, sd34

burcin commented 12 years ago

Reviewer: Burcin Erocal

malb commented 12 years ago

Changed reviewer from Burcin Erocal to Burcin Erocal, Simon King, Martin Albrecht

malb commented 12 years ago

comment:25

Burcin's patch looks good. Thus, giving it a positive review. I'm also running doctests again against 4.7.2.alpha3.

malb commented 12 years ago

Description changed:

--- 
+++ 
@@ -32,3 +32,4 @@
 * **Apply** [attachment: trac_4260-matrix-modn-docs.patch](https://github.com/sagemath/sage-prod/files/10642187/trac_4260-matrix-modn-docs.patch.gz)
 * **Apply** [attachment: trac_4260_more_doctests.patch](https://github.com/sagemath/sage-prod/files/10642188/trac_4260_more_doctests.patch.gz)
 * **Apply** [attachment: trac_4260_echelonformdomain.patch](https://github.com/sagemath/sage-prod/files/10642190/trac_4260_echelonformdomain.patch.gz)
+* **Apply** [attachment: trac_4260-minor_fixes.patch](https://github.com/sagemath/sage-prod/files/10642191/trac_4260-minor_fixes.patch.gz)

malb commented 12 years ago

comment:27

Doctests indeed pass on sage.math.

jdemeyer commented 12 years ago

Changed work issues from improve echelonize to none

jdemeyer commented 12 years ago

comment:29

The spkg could do with some cleanup:

What is the purpose of spkg-rebuild? If it is not used, remove it. If it is used, document it.
spkg-debian and the dist directory should be removed. They are leftovers for Debian, but these are now being removed from every spkg.
Why is "linbox" in .hgignore?
The file patches/commentator.C lacks a corresponding .patch file.

Optional:

Make spkg-install use patch for patching.

jdemeyer commented 12 years ago

Work Issues: cleanup spkg

malb commented 12 years ago

comment:30

What is the purpose of spkg-rebuild? If it is not used, remove it. If it is used, document it.

removed

spkg-debian and the dist directory should be removed. They are leftovers for Debian, but these are now being removed from every spkg.

removed

Why is "linbox" in .hgignore?

removed

The file patches/commentator.C lacks a corresponding .patch file.

added

Make spkg-install use patch for patching.

left for another time

malb commented 12 years ago

Description changed:

--- 
+++ 
@@ -26,7 +26,7 @@
 Time: 1.560

- Install http://sage.math.washington.edu/home/malb/spkgs/linbox-1.1.6.p4.spkg + Install http://sage.math.washington.edu/home/malb/spkgs/linbox-1.1.6.p5.spkg

Apply attachment: trac_4260-linbox_default.patch
Apply attachment: trac_4260-dense_ctypes_template.patch
Apply attachment: trac_4260-matrix-modn-docs.patch

jdemeyer commented 12 years ago

Changed work issues from cleanup spkg to none

jdemeyer commented 12 years ago

Merged: sage-4.7.3.alpha0

sagemath / sage

use LinBox as native matrix representation for dense matrices over GF(p) #4260