csEtree.js 1.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869
  1. /**
  2. * Computes the elimination tree of Matrix A (using triu(A)) or the
  3. * elimination tree of A'A without forming A'A.
  4. *
  5. * @param {Matrix} a The A Matrix
  6. * @param {boolean} ata A value of true the function computes the etree of A'A
  7. *
  8. * Reference: http://faculty.cse.tamu.edu/davis/publications.html
  9. */
  10. export function csEtree(a, ata) {
  11. // check inputs
  12. if (!a) {
  13. return null;
  14. }
  15. // a arrays
  16. var aindex = a._index;
  17. var aptr = a._ptr;
  18. var asize = a._size;
  19. // rows & columns
  20. var m = asize[0];
  21. var n = asize[1];
  22. // allocate result
  23. var parent = []; // (n)
  24. // allocate workspace
  25. var w = []; // (n + (ata ? m : 0))
  26. var ancestor = 0; // first n entries in w
  27. var prev = n; // last m entries (ata = true)
  28. var i, inext;
  29. // check we are calculating A'A
  30. if (ata) {
  31. // initialize workspace
  32. for (i = 0; i < m; i++) {
  33. w[prev + i] = -1;
  34. }
  35. }
  36. // loop columns
  37. for (var k = 0; k < n; k++) {
  38. // node k has no parent yet
  39. parent[k] = -1;
  40. // nor does k have an ancestor
  41. w[ancestor + k] = -1;
  42. // values in column k
  43. for (var p0 = aptr[k], p1 = aptr[k + 1], p = p0; p < p1; p++) {
  44. // row
  45. var r = aindex[p];
  46. // node
  47. i = ata ? w[prev + r] : r;
  48. // traverse from i to k
  49. for (; i !== -1 && i < k; i = inext) {
  50. // inext = ancestor of i
  51. inext = w[ancestor + i];
  52. // path compression
  53. w[ancestor + i] = k;
  54. // check no anc., parent is k
  55. if (inext === -1) {
  56. parent[i] = k;
  57. }
  58. }
  59. if (ata) {
  60. w[prev + r] = k;
  61. }
  62. }
  63. }
  64. return parent;
  65. }