# Reconstruction on large sparse array

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[[ ;; ]] ]];

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[[ ;; ]]]];

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!

• It's hard to say anything since we don't have your SparseArray[]. Could you consider using a matrix from e.g. ExampleData["Matrix"] as an example? Jul 22, 2022 at 16:32
• Hi, ExampleData[{"Matrix", "685BUS"}] looks similar to what I'm dealing with, only my matrix is much bigger. I also just fixed the numbering of the From vectors in my original post, since I falsely used From1 in all lines. Please pardon me! Jul 22, 2022 at 17:48
• Could you actually give a small example with specific From and To lists. I think I get what you're saying, but it'd be nice to have a test case to verify against. So, like provide a starting matrix, a "from" list, a "to" list, and then the expected result matrix. It doesn't seem like the logic will depend on the fact that the array is sparse. It also seems like the last bit about a "take" list is separate, but maybe more concreteness there would help. Jul 23, 2022 at 22:46
• Hi, I just added a 10x10 example in the original post (see above), I marked the corresponding part with todays date. The idea/procedure is not depending on sparseness, as you said, but I suspect that the solution for speeding up this procedure might be... Thanks a lot for dealing! Jul 25, 2022 at 8:59

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]},
+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.

• Hi, your approach does the correct thing with my Kmat. My large sparse array Kmat has 996336 non-zero entries. The AbsoluteTiming of the reconstruction is reduced from 0.12 sec to 0.05 sec, which is a nice speedup! I was hoping to get it reduced to a negligible computation time, but perhaps that's too much to ask for since I'm dealing wth huge data... Your first point in the disclaimers is correct: The Takes is the complement of all equation numbers and the Froms, the Tos stay in the new matrix. Thanks a lot, @user293787! Jul 27, 2022 at 13:27
• @Frida Thank you for the feedback. That is only a moderate speedup. I am assuming that in your case X is always the same and only Kmat keeps changing. Perhaps one could exploit that, do some preprocessing on X and so on. But it is probably only worth doing if it is a bottleneck in your computation. Jul 27, 2022 at 13:35