algorithm for index numbers of triangular matrix coefficients

Question

Hi all,

I think this must be simple but I can't get it right...

I have an MxM triangular matrix, the coefficients of which are stored in a vector, row by row. For example:

M =   [ m00 m01 m02 m03 ] 
      [     m11 m12 m12 ]
      [         m22 m23 ]
      [             m33 ]

is stored as

coef[ m00 m01 m02 m03 m11 m12 m13 m22 m23 m33 ]

Now I'm looking for a non-recursive algorithm that gives me for matrix size M and coefficient array index i

unsigned int row_index(i,M)

and

unsigned int column_index(i,M)

of the matrix element that it refers to. So, row_index(9,4) == 3, column_index(7,4) == 2 etc. if the index counting is zero-based.

EDIT: Several replies using an iteration have been given. Does anyone know of algebraic expressions?

Answer 1

Took me some time to understand what you needed! :)

unsigned int row_index(int i, int m)
{
    int iCurrentRow = 0;
    int iTotalItems = 0;
    for(int j = m; j > 0; j--)
    {
     iTotalItems += j;

     if( (i+1) <= iTotalItems)
      return iCurrentRow;

     iCurrentRow ++;
    }

    return -1; // Not checking if "i" can be in a MxM matrix.
}

Sorry forgot the other function.....

unsigned int column_index(int i, int m)
{
    int iTotalItems = 0;
    for(int j = m; j > 0; j--)
    {
     iTotalItems += j;

     if( (i+1) <= iTotalItems)
      return m - (iTotalItems - i);
    }

    return -1; // Not checking if "i" can be in a MxM matrix.
}

Answer 2

It has to be that

i == col + row*(M-1)-row*(row-1)/2

So one way to find col and row is to iterate over possible values of row:

for(row = 0; row < M; row++){
  col = i - row*(M-1)-row*(row-1)/2
  if (row <= col < M) return (row,column);
}

This is at least non-recursive, I don't know if you can do this without iteration.

As can be seen from this and other answers, you pretty much have to calculate row in order to calculate column, so it could be smart to do both in one function.

Answer 3

There may be a clever one liner for these, but (minus any error checking):

unsigned int row_index( unsigned int i, unsigned int M ){
    unsigned int row = 0;
    unsigned int delta = M - 1;
    for( unsigned int x = delta; x < i; x += delta-- ){
        row++;
    }
    return row;
}

unsigned int column_index( unsigned int i, unsigned int M ){
    unsigned int row = 0;
    unsigned int delta = M - 1;
    unsigned int x;
    for( x = delta; x < i; x += delta-- ){
        row++;
    }
    return M + i - x - 1;
}

Answer 4

Here's an algebraic (mostly) solution:

unsigned int row_index( unsigned int i, unsigned int M ){
    double m = M;
    double row = (-2*m - 1 + sqrt( (4*m*(m+1) - 8*(double)i - 7) )) / -2;
    if( row == (double)(int) row ) row -= 1;
    return (unsigned int) row;
}


unsigned int column_index( unsigned int i, unsigned int M ){
    unsigned int row = row_index( i, M);
    return  i - M * row + row*(row+1) / 2;
}

EDIT: fixed row_index()

Answer 5

Hi, I thought a little bit and I got the following result. Note that you get both row and columns on one shot.

Assumptions: Rows start a 0. Columns start at 0. Index start at 0

Notation

N = matrix size (was M in the original problem)

m = Index of the element

The psuedo code is

function ind2subTriu(m,N)
{
  d = 0;
  i = -1;
  while d < m
  {
    i = i + 1
    d = i*(N-1) - i*(i-1)/2
  }
  i0 = i-1;
  j0 = m - i0*(N-1) + i0*(i0-1)/2 + i0 + 1;
  return i0,j0
}

And some octave/matlab code

function [i0 j0]= ind2subTriu(m,N)
 I = 0:N-2;
 d = I*(N-1)-I.*(I-1)/2;
 i0 = I(find (d < m,1,'last'));
 j0 = m - d(i0+1) + i0 + 1;

What do you think?