Skip to content

Commit 92be1f9

Browse files
author
Release Manager
committed
gh-40773: Faster determinant for matrices over gf2e (M4RIE) This computes the determinant for matrices over `GF(2**e)` by directly relying on the PLE decomposition of M4RIE. This provides a significant speed-up over the generic implementation used until now for this class. Before: ``` sage: mat = matrix.random(GF(2**4), 200, 200) sage: %time dd = mat.det() CPU times: user 706 ms, sys: 2.26 ms, total: 708 ms Wall time: 705 ms sage: mat = matrix.random(GF(2**8), 200, 200) sage: %time dd = mat.det() CPU times: user 746 ms, sys: 890 µs, total: 747 ms Wall time: 744 ms sage: mat = matrix.random(GF(2**15), 200, 200) sage: %time dd = mat.det() CPU times: user 1.37 s, sys: 4.74 ms, total: 1.38 s Wall time: 1.37 s sage: mat = matrix.random(GF(2**4), 500, 500) sage: %time dd = mat.det() CPU times: user 10.8 s, sys: 23.8 ms, total: 10.8 s Wall time: 10.7 s sage: mat = matrix.random(GF(2**8), 500, 500) sage: %time dd = mat.det() CPU times: user 11.6 s, sys: 22.1 ms, total: 11.6 s Wall time: 11.6 s sage: mat = matrix.random(GF(2**15), 500, 500) sage: %time dd = mat.det() CPU times: user 21.7 s, sys: 37 ms, total: 21.7 s Wall time: 21.7 s ``` After: ``` sage: mat = matrix.random(GF(2**4), 200, 200) sage: %time dd = mat.det() CPU times: user 275 µs, sys: 68 µs, total: 343 µs Wall time: 356 µs sage: mat = matrix.random(GF(2**8), 200, 200) sage: %time dd = mat.det() CPU times: user 1.1 ms, sys: 21 µs, total: 1.12 ms Wall time: 1.13 ms sage: mat = matrix.random(GF(2**15), 200, 200) sage: %time dd = mat.det() CPU times: user 81.2 ms, sys: 8.98 ms, total: 90.1 ms Wall time: 89.7 ms sage: mat = matrix.random(GF(2**4), 500, 500) sage: %time dd = mat.det() CPU times: user 1.84 ms, sys: 0 ns, total: 1.84 ms Wall time: 1.86 ms sage: mat = matrix.random(GF(2**8), 500, 500) sage: %time dd = mat.det() CPU times: user 5.23 ms, sys: 39 µs, total: 5.27 ms Wall time: 5.44 ms sage: mat = matrix.random(GF(2**15), 500, 500) sage: %time dd = mat.det() CPU times: user 719 ms, sys: 17.6 ms, total: 737 ms Wall time: 732 ms ``` URL: #40773 Reported by: Vincent Neiger Reviewer(s): user202729, Vincent Neiger
2 parents d18ce74 + 8ee83ca commit 92be1f9

File tree

2 files changed

+86
-2
lines changed

2 files changed

+86
-2
lines changed

src/sage/libs/m4rie.pxd

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -5,7 +5,7 @@
55
# http://www.gnu.org/licenses/
66
##############################################################################
77

8-
from sage.libs.m4ri cimport mzd_t, m4ri_word
8+
from sage.libs.m4ri cimport mzd_t, m4ri_word, mzp_t
99

1010

1111
cdef extern from "m4rie/m4rie.h":
@@ -109,6 +109,8 @@ cdef extern from "m4rie/m4rie.h":
109109

110110
size_t mzed_echelonize_ple(mzed_t *, size_t)
111111

112+
int mzed_ple(mzed_t *A, mzp_t *P, mzp_t *Q)
113+
112114
#cdef extern from "m4rie/strassen.h":
113115
mzed_t *mzed_mul_strassen(mzed_t *, mzed_t *, mzed_t *, size_t cutoff)
114116

src/sage/matrix/matrix_gf2e_dense.pyx

Lines changed: 83 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -95,7 +95,7 @@ from sage.misc.randstate cimport randstate, current_randstate
9595
from sage.matrix.matrix_mod2_dense cimport Matrix_mod2_dense
9696
from sage.matrix.args cimport SparseEntry, MatrixArgs_init
9797

98-
from sage.libs.m4ri cimport m4ri_word, mzd_copy
98+
from sage.libs.m4ri cimport m4ri_word, mzd_copy, mzp_t, mzp_init, mzp_free
9999
from sage.libs.m4rie cimport *
100100
from sage.libs.m4rie cimport mzed_t
101101

@@ -1602,6 +1602,88 @@ cdef class Matrix_gf2e_dense(matrix_dense.Matrix_dense):
16021602
mzed_cling(self._entries, v)
16031603
mzd_slice_free(v)
16041604

1605+
def determinant(self):
1606+
"""
1607+
Return the determinant of this matrix.
1608+
1609+
Relies directly on M4RIE's PLE decomposition, and incidentally caches
1610+
the rank of ``self``.
1611+
1612+
EXAMPLES::
1613+
1614+
sage: gf4.<z> = GF(4)
1615+
sage: mat = matrix(gf4, 2, 2, [[z + 1, z + 1], [z, 1]])
1616+
sage: mat
1617+
[z + 1 z + 1]
1618+
[ z 1]
1619+
sage: mat.determinant()
1620+
z
1621+
sage: gf256.<t> = GF(2**8)
1622+
sage: mat = matrix(gf256, 3, 3, [[1, t, t**2],
1623+
....: [t**2, 1, t],
1624+
....: [t, t**2, 1]])
1625+
sage: mat.determinant()
1626+
t^6 + 1
1627+
sage: mat = matrix(gf256, 3, 3, [[1, t, t**2],
1628+
....: [t**2, 1, t],
1629+
....: [t**2 + 1, t + 1, t**2 + t]])
1630+
sage: mat.determinant()
1631+
0
1632+
1633+
Non-square matrices and the `0 \times 0` matrix are taken care of::
1634+
1635+
sage: matrix(gf4, 0, 0).determinant()
1636+
1
1637+
sage: matrix(gf4, 3, 2).determinant()
1638+
Traceback (most recent call last):
1639+
...
1640+
ValueError: self must be a square matrix
1641+
"""
1642+
cdef size_t m = self._nrows
1643+
1644+
if m != self._ncols:
1645+
raise ValueError("self must be a square matrix")
1646+
if m == 0:
1647+
return self._one
1648+
1649+
x = self.fetch('det')
1650+
if x is not None:
1651+
return x
1652+
1653+
cdef mzed_t * A = mzed_copy(NULL, self._entries)
1654+
cdef mzp_t * P = mzp_init(m)
1655+
cdef mzp_t * Q = mzp_init(m)
1656+
1657+
sig_on()
1658+
cdef int r = mzed_ple(A, P, Q)
1659+
sig_off()
1660+
1661+
self.cache('rank', r)
1662+
1663+
if r < m:
1664+
mzp_free(P)
1665+
mzp_free(Q)
1666+
mzed_free(A)
1667+
self.cache('det', self._zero)
1668+
return self._zero
1669+
1670+
cdef Cache_base cache = <Cache_base> self._base_ring._cache
1671+
1672+
# characteristic 2, so det(P) == det(Q) == 1
1673+
cdef Py_ssize_t i
1674+
cdef int elt
1675+
cdef det = self._one
1676+
for i from 0 <= i < m:
1677+
elt = mzed_read_elem(A, i, i)
1678+
det = det * cache.fetch_int(elt)
1679+
1680+
mzp_free(P)
1681+
mzp_free(Q)
1682+
mzed_free(A)
1683+
1684+
self.cache('det', det)
1685+
return det
1686+
16051687

16061688
def unpickle_matrix_gf2e_dense_v0(Matrix_mod2_dense a, base_ring, nrows, ncols):
16071689
r"""

0 commit comments

Comments
 (0)