Reconstruction on large sparse array

Question

I'm struggling with a reconstruction of a large sparse array. What I'm trying to do is to rearrange a 2 dimensional array (respresenting a system matrix for a given system of equations) such that I add up certain columns and rows and only take a submatrix of the original one.

To be more specific:

Let Kmat be a large sparse array with the dimensions 44000x44000. The vectors From1, From2 and From3 give specific row/column numbers from where I want to take the values, which are added to the row/column numbers given in the vectors To1, To2 and To3. The dimensions of these vectors are much smaller, e.g. 300 to 400 entries.

Therefore, I want to perform:

(* adding rows *)    
Kmat[[;; , To1[[ ;; ]] ]] += Kmat[[;; , From1[[ ;; ]] ]];
Kmat[[;; , To2[[ ;; ]] ]] += Kmat[[;; , From2[[ ;; ]] ]];
Kmat[[;; , To3[[ ;; ]] ]] += Kmat[[;; , From3[[ ;; ]] ]];

(* adding columns *)
Kmat[[To1[[ ;; ]] ]] += Kmat[[From1[[ ;; ]] ]];
Kmat[[To2[[ ;; ]] ]] += Kmat[[From2[[ ;; ]] ]];
Kmat[[To3[[ ;; ]] ]] += Kmat[[From3[[ ;; ]] ]];

And then I only want to take a submatrix defined by the indices in the vector take with the dimensions around 43000:

Kmatnew = Kmat[[takes,takes]];

This all takes only about 0.1 seconds on my machine, which I would consider to be not expensive at all, but I repeat this a few thousand times which definitely adds up in my computation... I aim to do very efficient computing and that is definitely one of the bottlenecks of my code.

My first idea:

Try to use Map in some way, but I was not succesful in findiing a solution using Map...

My second idea:

Define a compile function with Compile, but unfortunately that does not work with sparse arrays.

I'd be very happy about some useful hints and ideas...

########################################

Edit 25th July 2022: Add a small example to illustrate:

Let's take an example in the dimensions 10x10:

dim = 10;
Kmat = Table[RandomInteger[10], {i, dim}, {j, dim}];

This results in e.g.

   Kmat = 
{{10, 10, 1, 2, 4, 3, 8, 2, 6, 6}, 
{7, 0, 5, 0, 7, 2, 3, 6, 6, 3}, 
{1,4, 3, 9, 8, 4, 8, 2, 5, 2}, 
{0, 5, 0, 5, 7, 3, 0, 8, 10, 0}, 
{10, 6, 1, 6, 7, 1, 5, 6, 0, 1}, 
{2, 9, 7, 3, 5, 7, 6, 8, 4, 9}, 
{5, 4, 5, 0, 5, 1, 1, 2, 2, 10}, 
{2, 9, 0, 3, 10, 1, 9, 10, 0, 8}, 
{0, 7, 4, 10, 5, 7, 8, 0, 9, 10}, 
{10, 10, 3, 9, 1, 1, 9, 3, 0, 4}};

We want to add the rows/columns in From to the rows/colums in To. Then, we take whats left if we leave out the Froms.

To1 = {1, 3}; To2 = {2, 4}; To3 = {1};
From1 = {6, 8}; From2 = {7, 9}; From3 = {7};
Takes = Complement[Range[dim], Union[From1, From2, From3]];

Now we perform the upper procedure to add:

(* adding rows *)
Kmat[[;; , To1[[ ;; ]]]] += Kmat[[;; , From1[[ ;; ]]]];
Kmat[[;; , To2[[ ;; ]]]] += Kmat[[;; , From2[[ ;; ]]]];
Kmat[[;; , To3[[ ;; ]]]] += Kmat[[;; , From3[[ ;; ]]]];

(* adding columns *)
Kmat[[To1[[ ;; ]]]] += Kmat[[From1[[ ;; ]]]];
Kmat[[To2[[ ;; ]]]] += Kmat[[From2[[ ;; ]]]];
Kmat[[To3[[ ;; ]]]] += Kmat[[From3[[ ;; ]]]];

And this should result in

    Kmat = 
{{43, 38, 25, 17, 14, 11, 15, 12, 12, 25}, 
{19, 8, 18, 8, 12, 3, 4, 8, 8, 13}, 
{25, 30, 15, 17, 18, 5, 17, 12, 5, 10}, 
{18, 20, 12, 34, 12, 10, 8, 8, 19, 10}, 
{16, 11, 7, 6, 7, 1, 5, 6, 0, 1},
{15, 15, 15, 7, 5, 7, 6, 8, 4, 9}, 
{7, 5, 7, 2, 5, 1, 1, 2, 2, 10}, 
{12, 18, 10, 3, 10, 1, 9, 10, 0, 8}, 
{15, 15, 4, 19, 5, 7, 8, 0, 9, 10}, 
{20, 19, 6, 9, 1, 1, 9, 3, 0, 4}}

Then, we take the Takes:

Kmatnew = Kmat[[Takes, Takes]];

and Kmatnew results in

Kmatnew = 
{{43, 38, 25, 17, 14, 25}, 
{19, 8, 18, 8, 12, 13}, 
{25, 30, 15, 17, 18, 10}, 
{18, 20, 12, 34, 12, 10}, 
{16, 11, 7, 6, 7, 1}, 
{20, 19, 6, 9, 1, 4}}

Sorry for the long post!

Answer 1

Disclaimers. This answer comes with two disclaimers:

• I checked correctness only for the dim=10 toy example given by OP. In the toy example, Takes is the complement of all From, and the To are contained in Takes, and so on. I do not know if that is supposed to be typical. Also, I have not clarified what exact assumptions one needs for this answer to give the correct result.
• I have not tested speed. Whether my approach is useful for OP in terms of speed could depend, for example, on how many nonzero the input matrix Kmat typically has, which OP has not specified. Without a typical example at hand, it did not seem a good idea to test speed.

Answer. Given dim, From1, From2, From3, To1, To2, To3 and Takes we can define the following auxiliary sparse matrix

X = With[{To=Join[To1,To2,To3],From=Join[From1,From2,From3]},
        (SparseArray[MapThread[(List[##]->1)&,{To,From}],{dim,dim}]
           +IdentityMatrix[dim,SparseArray])[[Takes,;;]]
];

It should have on the order of dim many nonzero entries.

Then, given Kmat, we calculate the new matrix using two sparse matrix multiplications:

Kmatnew = X.Kmat.Transpose[X]

Example. As already mentioned, correctness was only checked in the dim=10 toy example, and only for the particular From1, From2, From3, To1, To2, To3 and Takes given by OP. In this case, the auxiliary matrix X is the sparse array version of

{{1,0,0,0,0,1,1,0,0,0},
 {0,1,0,0,0,0,1,0,0,0},
 {0,0,1,0,0,0,0,1,0,0},
 {0,0,0,1,0,0,0,0,1,0},
 {0,0,0,0,1,0,0,0,0,0},
 {0,0,0,0,0,0,0,0,0,1}}

In this case the result is correct for all $10\times 10$ matrices Kmat.