Apologies if this isn't the appropriate forum for this question.

I have a set of elements that I need to iterate over as part of a modeling workflow. The elements exists over a set of dimensions (i, j, k, l). Due to model constraints, most of the elements within the set are not considered (value = 0). Therefore, looping over the set in the straightforward way is inefficient. Eg:

#: Not a good way to do it:
for i in range(n_i):
    for j in range(n_j):
        for k in range(n_k):
            for l in range(n_l):

Since most of the elements are 0, I should be able to eliminate them from the inner loops given the index values from the outer loops. As an example, let's say that n_j = 10, but for i=1, only j = 2,3 are considered by the model. Then for i=1, I should only have to iterate over j=2,3 which avoids ~80% of the iterations on that loop.

What I would like to do is rewrite my iterable range as a function of the outer loop parameters. Something like:

N_i = I(j,k,l)
N_j = J(i,k,l)
N_k = K(i,j,l)
N_l = L(i,j,k)

#: Want to do something like this
for i in I():                        #: loop over all i
    for j in J(i):                   #: loop over j for i = #
        for k in K(i,j):             #: loop over k for i = #, j = #
            for l in L(i,j,k):       #: loop over l for i = #, j = #, k = #

where N_i ... N_l are the subset of that iterable set that is != 0, given one or more of the i,j,k,l.

I could imaging building a nested dictionary to look it up, but I think i would need a different dictionary for every possible combination of loop order i,j,k,l -> i,l,j,k -> l,i,j,k etc...

My question is, what's an efficient way to structure the data & write these functions in order to achieve this end? Also if it matters, I will iterate over the dimensions in a different order for different aspects of the model.

I'm doing this all in Python so answers that address the python implementation would be great.

  • $\begingroup$ Even with 2-dimensional sparse matrices, you can typically iterate efficiently in one loop order (i.e. Row Major i->j), otherwise you're searching. If that's not acceptable, then maybe you're better off storing 4! = 24 copies of the tensor! $\endgroup$ – Charlie S Aug 14 '20 at 20:09

Start with how this would look like for a sparse matrix, i.e., only two indices. In that case, data is generally stored in "Compressed Row Storage" (CSR) or "Compressed Column Storage" (CSC). If you understand how data is stored in these formats, you will also understand how to store data for the higher-dimensional cases -- these are generally called "sparse tensors" and you will find a substantial amount of literature on that case.

  • $\begingroup$ The data all fits in memory, so I'm not particularly worried about storage efficiency. I just need to be able to loop through it quickly. I'll check out sparse tensors and thanks for the lead. $\endgroup$ – Sledge Aug 14 '20 at 12:43
  • 1
    $\begingroup$ It's not about storage efficiency alone. If you already store where the nonzeros are, you might as well only store the nonzeros. $\endgroup$ – Wolfgang Bangerth Aug 14 '20 at 13:33
  • $\begingroup$ I've got that part figured out. What I still struggle with is how to iterating over my data without looping over every i,j,k,l. $\endgroup$ – Sledge Aug 14 '20 at 15:05
  • 1
    $\begingroup$ Think about how you would do that for CSR: You'd loop over all $i$, and for each $i$ you've already go stored which indices $j$ exist (i.e., which entries $A_{ij}\neq 0$). $\endgroup$ – Wolfgang Bangerth Aug 14 '20 at 18:23

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.