Add faster implementations of matrix_from_* constructors to GF(2) and GF(2^e)

[Biffo89]

Dense matrices over GF(2) and GF(2^e) currently have slow performance for the following methods:

• matrix_from_rows()
• matrix_from_columns()
• matrix_from_rows_and_columns()

They are using the default implementation in matrix1.pyx which involves conversion of each element every time it is read or written. This PR puts direct naive implementations for these methods into both classes.

📝 Checklist

• The title is concise and informative.
• The description explains in detail what this PR is about.
• I have linked a relevant issue or discussion.
• I have created tests covering the changes.
• I have updated the documentation and checked the documentation preview.

[user202729]

while this should be correct (haven't read in detail), wouldn't it be better to consolidate the fast implementation into shared template somehow and each subclass only need to override one or two methods?

Since virtual method are supported in Python, we already have get_unsafe and set_unsafe, we could also provide a virtual method A.copy_from_unsafe(i, j, B, k, l) that can be called generically.

The cost of dereferencing a single function pointer ought to be small enough, and the function pointer is known to be the same over all iterations.

@tscrim any opinion?

[Biffo89]

How's this for a replacement?

In matrix_generic_dense.pyx:

    cdef copy_from_unsafe(self, Py_ssize_t iDst, Py_ssize_t jDst, src, Py_ssize_t iSrc, Py_ssize_t jSrc):
        self._entries[iDst*self._ncols + jDst] = src._entries[iSrc*self._ncols + jSrc]

    def transpose(self):
        (nc, nr) = (self.ncols(), self.nrows())
        cdef Matrix_dense trans
        trans = self.new_matrix(nrows = nc, ncols = nr,
                                copy=False, coerce=False)

        cdef Py_ssize_t i, j
        for j from 0<= j < nc:
            for i from 0<= i < nr:
                trans.copy_from_unsafe(j, i, self, i, j)

        if self._subdivisions is not None:
            row_divs, col_divs = self.subdivisions()
            trans.subdivide(col_divs, row_divs)
        return trans

In matrix1.pyx:

    def matrix_from_columns(self, columns):
        columns = PySequence_Fast(columns, "columns is not iterable")
        cdef Py_ssize_t ncols = len(columns)

        # Construct new matrix
        cdef Matrix A = self.new_matrix(ncols=ncols)
        cdef Py_ssize_t i, j, col
        for j, col in enumerate(columns):
            if col < 0 or col >= self._ncols:
                raise IndexError("column index out of range")
            for i in range(self._nrows):
                A.copy_from_unsafe(i, j, self, i, col)
        return A

    def matrix_from_rows(self, rows):
        rows = PySequence_Fast(rows, "rows is not iterable")
        cdef Py_ssize_t nrows = len(rows)

        # Construct new matrix
        cdef Matrix A = self.new_matrix(nrows=nrows)
        cdef Py_ssize_t i, j, row
        for i, row in enumerate(rows):
            if row < 0 or row >= self._nrows:
                raise IndexError("row index out of range")
            for j in range(self._ncols):
                A.copy_from_unsafe(i, j, self, row, j)
        return A

    def matrix_from_rows_and_columns(self, rows, columns):
        rows = PySequence_Fast(rows, "rows is not iterable")
        columns = PySequence_Fast(columns, "columns is not iterable")

        cdef Py_ssize_t ncols = len(columns)
        cdef Py_ssize_t nrows = len(rows)

        # Check whether column indices are valid
        cdef Py_ssize_t i, j, row, col
        for col in columns:
            if col < 0 or col >= self._ncols:
                raise IndexError("column index out of range")

        # Construct new matrix
        cdef Matrix A = self.new_matrix(nrows=nrows, ncols=ncols)
        for i, row in enumerate(rows):
            if row < 0 or row >= self._nrows:
                raise IndexError("row index out of range")
            for j, col in enumerate(columns):
                A.copy_from_unsafe(i, j, self, row, col)
        return A

In matrix_gf2e_dense.pyx (for example):

    cdef copy_from_unsafe(self, Py_ssize_t iDst, Py_ssize_t jDst, src, Py_ssize_t iSrc, Py_ssize_t jDst):
        mzed_write_elem(self._entries, iDst, jDst, mzed_read_elem(src, iSrc, jSrc))

[xcaruso]

I just have two comments:

• in the generic implementation, you can use get_unsafe and set_unsafe,
• in the doctest of copy_from_unsafe, you should insist on the fact that the method assumes that src has also type MatrixGF2e_dense.

[Biffo89]

I'm working on implementing copy_from_unsafe now. Given the #40423 and #40425 can both use this, would it be worth updating those with copy_from_unsafe and making the two PRs dependent on this one?

[xcaruso]

Yes, sure!

[user202729]

The poblem with {get,set}_unsafe is that for it to work in generic code, the return value need to be converted to a Python object, which is not entirely efficient in many cases (e.g. if it can be kept as a C type). copy_from_unsafe doesn't suffer from this problem.

[Biffo89]

The poblem with {get,set}_unsafe is that for it to work in generic code, the return value need to be converted to a Python object, which is not entirely efficient in many cases (e.g. if it can be kept as a C type). copy_from_unsafe doesn't suffer from this problem.

Yes, a generic inefficient copy_from_unsafe can be used as a default implementation, and in practice the vast majority of classes can override this.

[Biffo89]

I've added a copy_from_unsafe to all files with a get_unsafe/set_unsafe implementation except matrix_window.pyx which I wasn't sure how to handle. All tests are passing.

#40425 should be obsolete with this PR, #40423 is still necessary as it is a bug fix.

[xcaruso]

Thanks! Everything looks good to me.