Performance issue with sparse matrix method "TransposeThisAndMultiply"


(richard reader) #1

Hello MathNet

This is my first time posting on this site, so please let me know if I miss anything that I should add to this message, or if the information I give is not the way it should be.

I have noticed that the method “Matrix.TransposeThisAndMultiply(Vector rightSide)” is extremely slow for sparse matrices. This is due to the fact that the sparse matrix implementations do not override the the base class method (i.e. Matrix.DoTransposeThisAndMultiply(Vector rightSide)). This base class method calls the function “At(i, j)” which searches the sparse matrices for elements at a given (row, column) position - this is a very slow operation for sparse matrices.

To avoid this problem, users could not use this function at all, but instead transpose their sparse matrices, and then call “Multiply”. Alternatively, the SparseMatrix classes could override TransposeThisAndMultiply(Vector rightSide) with an implementation suitable for sparse matrices stored in CSR format. I have pasted below a possible new implementation.

For a sparse matrix of 20,000 rows, 5,000 columns and 100 non-zero values per row, I get the following timings:

  1. using base class implementation of DoTransposeThisAndMultiply(Vector rightSide) takes around 15 seconds to multiply a single vector by this matrix

  2. transposing the matrix, then calling Multiply, takes around 0.15 seconds, but uses more memory as the transposed matrix must be created and stored.

  3. the suggested implementation below takes around 0.17 seconds, and has no additional memory requirements.

I would suggest the best route would be to incorporate the implementation below (or similar). Please let me know if the MathNet team/users agree, or if there is a better course of action to follow?

Thanks
Richard

Alternative Implementation which could be added to SparseMatrix classes:

    protected override void DoTransposeThisAndMultiply(Vector<double> rightSide, Vector<double> result)
    {
        var rowPointers = _storage.RowPointers;
        var columnIndices = _storage.ColumnIndices;
        var values = _storage.Values;

        for (var row = 0; row < RowCount; row++)
        {
            var startIndex = rowPointers[row];
            var endIndex = rowPointers[row + 1];

            if (startIndex == endIndex)
            {
                continue;
            }

            var rightSideValue = rightSide[row];
            for (var index = startIndex; index < endIndex; index++)
            {
                result[columnIndices[index]] += values[index] * rightSideValue;
            }
        }
    }