More sparse functionality.

[mmelnich]

I need C := op(A) + β·B for some operator in RandLAPACK. MKL has it.

[rileyjmurray]

Quick fix if our COO--> CSR/CSC matrix conversion handles duplicate (val, row, col) entries by summing vals.

• Define (vals, rows, cols) triplet of std::vectors.
• Iterate over nonzeros in op(A); adding them to (vals, rows, cols) as appropriate.
• Iterate over nonzeros in B; add beta*(the values) to (vals, rows, cols) as appropriate.
• Construct a COOMatrix that relies on vals.data(), rows.data(), cols.data()
• Construct a CSRMatrix (or CSCMatrix) from the COOMatrix.
• Return the compressed matrix.

If our conversion functions don't handle duplicate entries in COO then I advise we adjust them so they do.