Increasing the speed of sparse array

Question

I am trying to find a sparse matrix of the size Col x Rows = 16777216 x 1024 where all the elements are 0 except the ones given by col2 and row2. This is related to a previous question I asked: Link to the question

The answer by kglr solved my problem for smaller size of the sparse matrix but when I use bigger sizes of Rows, Col, col2, row2, as in the code given below, the computation goes on forever with no way of knowing when this would complete.

Col = Flatten[Outer[{#4, #3, #2, #1} &, Range[0, 1023, 1], Range[0, 31, 1], {1}, DeleteCases[Range[256, -256, -1], 0, Infinity], 1], 3];

Rows = Flatten[Outer[{#1, #2} &, DeleteCases[Range[-16, 16, 1], 0, Infinity], DeleteCases[Range[-16, 16, 1], 0, Infinity], 1], 1];

col2 = Flatten[Outer[{#4, #3, #2, #1} &, RandomChoice[Range[0, 1023, 1], 64], RandomChoice[Range[0, 31, 1], 16], {1}, RandomChoice[DeleteCases[Range[256, -256, -1], 0, Infinity], 61], 1], 3];

row2 = Rows;

positions = Tuples[{Flatten@Position[Col, Alternatives @@ col2], 
                 Flatten@Position[Rows, Alternatives @@ row2]}];

SparseMat = SparseArray[positions -> (f[Flatten[{Col[[#[[1]]]], Rows[[#[[2]]]]}]] & /@ 
                        positions), Length /@ {Col, Rows}];

My question is how can I speed up this computation?. Is there a faster way to do this?.

Answer 1

This uses packed arrays throughout and replaces Position with the significantly faster Nearest. (Position[Col, Alternatives @@ col2] needs about Length[Col] Length[col2] operations; Nearest[Col -> Automatic, col2, {1, 1}] needs only Log[Length[Col]] Length[col2] operations. So that's $10^{12}$ vs. $10^{6}$ operations.)

Moreover, compiling the function f such that it is threaded and parallelized over its first argument speeds up the evaluation f.

Admittedly, there are much faster ways to do this for the f provided by OP. But this should also accelerate more general functions.

First the compiled function

cf = Compile[{{X, _Integer, 1}, {Y, _Integer, 2}},
   Table[
    Plus[
     Compile`GetElement[X, 1],
     Compile`GetElement[X, 2],
     Compile`GetElement[X, 3],
     Compile`GetElement[X, 4],
     Compile`GetElement[Y, i, 1],
     Compile`GetElement[Y, i, 2]
     ], {i, 1, Length[Y]}],
   CompilationTarget -> "C",
   RuntimeAttributes -> {Listable},
   Parallelization -> True,
   RuntimeOptions -> "Speed"
   ];

Now, the restructured code runs through in about 6 seconds:

n = 1024;
Col = Tuples[{
    Range[0, n - 1, 1],
    Range[0, 31, 1],
    {1},
    DeleteCases[Range[256, -256, -1], 0, 1]
    }][[All, {4, 3, 2, 1}]];
Rows = Tuples[DeleteCases[Range[-16, 16, 1], 0, 1], 2];
nf = Nearest[Col -> Automatic];


col2 = Tuples[{
    RandomChoice[Range[0, n - 1, 1], 64],
    RandomChoice[Range[0, 31, 1], 16],
    {1},
    RandomChoice[DeleteCases[Range[256, -256, -1], 0, 1], 61]
    }][[All, {4, 3, 2, 1}]];
p1 = Flatten[Flatten[nf[col2]]];
positions = Tuples[{p1, Range[Length[Rows]]}];
values = Flatten[cf[Col[[p1]], Rows]];
SparseMat = SparseArray[positions -> values, Length /@ {Col, Rows}, 0];