1

I have the upper triangular portion of a matrix, with the main diagonal stored as a linear array, how can the (i,j) indices of a matrix element be extracted from the linear index of the array?

For example the linear array :[a0, a1, a2, a3, a4, a5, a6, a7, a8, a9, a10]is storage for the matrix

a0  a1  a2  a3
0   a4  a5  a6
0   0   a7  a8
0   0   0   a10

I've found solutions for this problem but without the main diagonal which is:

index = (n*(n-1)/2) - (n-i)*((n-i)-1)/2 + j - i - 1

And solution for the same problem but for a lower triangular matrix with the diagonal:

index = ((i + 1) * i / 2 + i).

Regards,

Bako
  • 313
  • 1
  • 15
  • Sounds like a homework question. What have you tried? What results did it produce? Can you modify the one without the main diagonal to work with the main diagonal? – mkasberg Sep 30 '16 at 14:37
  • I've tried this for without the diagonal : k = (n*(n-1)/2) - (n-i)*((n-i)-1)/2 + j - i - 1 And this for the lower triangualar with the diagonal :((I + 1) * I / 2 + J). – Bako Sep 30 '16 at 14:38
  • @Bako can you include this in the question along with code you have used. What you have provided here is the logic and not everyone looks into the comments. So, please edit the question. – The Apache Sep 30 '16 at 14:41
  • Do you have a mistake in “with the main diagonal stored as a linear array”? It appears that the linear *vector* stores the upper-triangular section, including the main diagonal? – Ahmed Fasih Sep 30 '16 at 15:54
  • @AhmedFasih No, I meant I store the upper triangular portion of a matrix including the main diagonal. – Bako Sep 30 '16 at 16:02
  • Got it! There’s a comma missing, but now I see – Ahmed Fasih Sep 30 '16 at 16:03
  • I can’t imagine how someone might think this is a homework question… a professor would have to be truly satanic to stoop to this level of pointlessness. Everyone who works with sparse, banded, or triangular arrays has to figure out this kind of thing if their language/library doesn’t… – Ahmed Fasih Sep 30 '16 at 16:46
  • @AhmedFasih He is right, It is actually a homework from my OOP class. I just didn't how to define this method, but I figured out in the end :D – Bako Oct 01 '16 at 14:49

3 Answers3

1

My solution might be equivalent to your’s, I haven’t checked:

index = N * i - ((i - 1) * i) / 2 + (j - i)

Here’s a complete Python test for it. I used Python because Numpy has triu_indices, which gives the upper-triangular indexes.

import numpy as np

def mksquare(N):
    """Make a square N by N matrix containing 0 .. N*N-1"""
    return np.arange(N * N).reshape(N, N)

def mkinds(N):
    """Return all triu indexes for N by N matrix"""
    return [(i,j) for i in range(N) for j in range(N) if i <= j]

def ij2linear(i, j, N):
    """Convert (i,j) 2D index to linear triu index for N by N array"""
    return N * i - ((i - 1) * i) // 2 + (j - i)

def test(N):
    """Make sure my `mkinds` works for given N"""
    arr = mksquare(N)
    vec = arr[np.triu_indices(N)]

    inds = mkinds(N)
    expected = [arr[i, j] for (i, j) in inds]

    actual = [vec[ij2linear(i, j, N)] for (i, j) in inds]

    return np.all(np.equal(actual, expected))

"""Run `test` for a bunch of `N`s and make sure they're all right"""
print(all(map(test, range(2, 20))))
# prints True

Worth a blog post explaining how to arrive at this conclusion, but this’ll do for now .

Community
  • 1
  • 1
Ahmed Fasih
  • 6,458
  • 7
  • 54
  • 95
0

I figured out the answer! It is :

index = (n*(n+1)/2) - (n-i)*((n-i)+1)/2 + j - i
Bako
  • 313
  • 1
  • 15
0

Triangular matrices generalize to any dimension. Suppose, for example, A[a,b,c,d,e] is non-zero only if 0<=a<=b<=c<=d<=e< N. We can compress A to a linear array X where A[a,b,c,d,e] = X[a+B[b]+C[c]+D[d]+E[e]]

Where B[b] = {b+1 choose 2}, C[c] = {c+2 choose 3}, D[d] = {d+3 choose 4} and E[e] = {e+4 choose 5}.

These "offset" arrays can be computed without using multiplication or division as follows:

B[0] = C[0] = D[0] = E[0] = 0;

for(int t = 1; t < N; t++)
 {
  B[t] = B[t-1]+t;
  C[t] = C[t-1]+B[t];
  D[t] = D[t-1]+C[t];
  E[t] = E[t-1]+D[t];
 }
Lawrence Larmore
  • 1