Intrusive access to matrix formats

[upsj]

Currently, I see no way how to modify the internal data structures of our matrix formats. By this I don't mean their values, but e.g. resizing the value array, as required for the SpGeMM and other methods modifying the sparsity pattern/nnz. As far as I can see, this is currently not possible from the public interface of the sparse formats? How would you extend the interface? Aside from resizing, the current pointer-interface should handle all use cases, so adding a combined resize_... (or set_nnz?) method for the internal arrays that can change independent of the dimensions seems like a simple solution. Other ideas?

[thoasm]

Do you really need that for the SpGeMM? This would happen from within the Csr matrix, so you would have access to all necessary attributes, and could also change them at will. Currently, it is also possible for the user to do a resize, albeit it is slightly more work than just calling resize: Creating a new matrix with the desired size, copying the values and then using move_to or copy_from to get it to the appropriate matrix.

The problem I see with a user accessable resize is that the matrix would be in an invalid state since you either remove non-zeros, resulting in row_ptrs pointing to elements which no longer exist (and the uncertainty which elements you want to remove), or adding values, leading to uncertainty in which column they should be put.

[hartwiganzt]

Without going into the discussion: SpMM = Sparse matrix Dense matrix SpGeMM = Sparse matrix Sparse matrix

[thoasm]

Good to know, thanks! I am pretty sure we were both talking about SpGeMM.

[upsj]

@thoasm with a free spmm kernel, this would not be possible, and I'd like to avoid having to separate it into two steps (symbolic SpGeMM + filling in the non-zeros) as not every implementation needs to have these two stages. My aim with this approach was to use as little temporary memory as possible by reusing the row_ptrs of the resulting matrix if they are of the correct size instead of allocating new storage and swapping the data array afterwards. But maybe this is a bit too over-zealous in terms of premature optimization?

Another idea would be to have a friend accessor class (MatrixBuilder?) that you need to use to intrusively modify these data structures? For example in each ParILUT iteration, we need to do at least three of these operations (SpGeMM, removing candidates for L and U)

[upsj]

Unless we start implementing our own very specialized sparsity-manipulating kernels, e.g. for SpGeMM with heavy focus on minimal memory usage, we can probably continue with the current two-step approach (counting nnz per row + prefix sum -> filling non-zeros per row) without this.