envelope.m

   1 function envelope = envelope(A)
   2   ## Compute the envelope of the matrix ``A``. The envelope of a matrix
   3   ## is defined as the set of indices,
   4   ##
   5   ##   E = { (i,j) : i < j, A(k,j) != 0 for some k <= i }
   6   ##
   7   if (!issymmetric(A) && !is_upper_triangular(A))
   8     ## The envelope of a matrix is only defined for U-T or symmetric
   9     ## matrices.
  10     envelope = {NA};
  11     return;
  12   end
  13
  14   ## Start with an empty result, and append to it as we find
  15   ## satisfactory indices.
  16   envelope = {};
  17
  18   for j = [ 1 : columns(A) ]
  19     ## Everything below the first non-zero element in a column will be
  20     ## part of the envelope. Since we're moving from top to bottom, we
  21     ## can simply set a flag indicating that we've found the first
  22     ## non-zero element. Thereafter, everything we encounter should be
  23     ## added to the envelope.
  24     found_nonzero = false;
  25
  26     for i = [ 1 : j-1 ]
  27       if (A(i,j) != 0)
  28         found_nonzero = true;
  29       end
  30
  31       if (found_nonzero)
  32         envelope{end+1} = [i,j];
  33       end
  34     end
  35
  36   end
  37
  38 end