How to speed up working with Sparse Arrays

Question

I make a SparseArray and then modify it by adding rows together. This seems to take a long time. Can it be speeded up?

na = 1000;
nb = 60000;
nc = 50;
c = 0.8;
data = Flatten[RandomReal[{-1, 1}, {na, na}]];
loc = RandomInteger[{1, 5386/2}, {na^2, 2}];
mat = SparseArray[loc -> data, {nb, nb}, 0.0];
ri = RandomInteger[{1, nb/2}, nc];

Now to do the adding of rows. This is what takes the time.

Timing[Do[mat[[n]] = mat[[n]] + c mat[[nb/2 + n]],
   {n, ri}];]

This takes 1.09 seconds on my machine Version; 11.3. Could it be done faster? Is there a workaround?

Thanks

Answer 1

Put as many operations as possible into one call:

mat[[ri]] += c mat[[nb/2 + ri]]; // AbsoluteTiming // First

0.010061

Edit

I also thought about a more general situation where rows enumerated by a list readidx are read from an input matrix, multiplied by a factor, and then added to rows enumerated by a writeidx. If readidx or writeidx contains any duplicates, one can create a lower triangular matrix to achieve the goal much quicker than with the Do loop:

SeedRandom[666];
nnz = 100000;
nb = 60000;
nc = 5000;

data = Flatten[RandomReal[{-1, 1}, nnz]];
loc = RandomInteger[{1, nb}, {nnz, 2}];
A = mat0 = SparseArray[loc -> data, {nb, nb}, 0.0];
readidx = RandomInteger[{1, nb}, nc];
writeidx = RandomInteger[{1, nb}, nc];
factors = RandomReal[{-1, 1}, nc];
DuplicateFreeQ[readidx]
DuplicateFreeQ[writeidx]

False

False

B1 = A;
Do[B1[[writeidx[[i]]]] += factors[[i]] A[[readidx[[i]]]], {i, 1, Length[readidx]}]; // AbsoluteTiming // First

B2 = With[{L = With[{spopt = SystemOptions["SparseArrayOptions"]},
        Internal`WithLocalSettings[
         SetSystemOptions["SparseArrayOptions" -> {"TreatRepeatedEntries" -> Total}],
         SparseArray[
          Transpose[{writeidx, readidx}] -> factors, 
          {Length[A], Length[A]}
         ],
         SetSystemOptions[spopt]]
        ]},
     A + L.A
     ]; // AbsoluteTiming // First
B1 == B2

18.1406

0.008159

True

Note that we have to SetSystemOptions["SparseArrayOptions" -> {"TreatRepeatedEntries" -> Total}] so that the matrix L gets assempled correctly, even if Transpose[{writeidx, readidx}] contains duplicates.