csCounts.js 3.4 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125
  1. import { factory } from '../../../utils/factory.js';
  2. import { csLeaf } from './csLeaf.js';
  3. var name = 'csCounts';
  4. var dependencies = ['transpose'];
  5. export var createCsCounts = /* #__PURE__ */factory(name, dependencies, _ref => {
  6. var {
  7. transpose
  8. } = _ref;
  9. /**
  10. * Computes the column counts using the upper triangular part of A.
  11. * It transposes A internally, none of the input parameters are modified.
  12. *
  13. * @param {Matrix} a The sparse matrix A
  14. *
  15. * @param {Matrix} ata Count the columns of A'A instead
  16. *
  17. * @return An array of size n of the column counts or null on error
  18. *
  19. * Reference: http://faculty.cse.tamu.edu/davis/publications.html
  20. */
  21. return function (a, parent, post, ata) {
  22. // check inputs
  23. if (!a || !parent || !post) {
  24. return null;
  25. }
  26. // a matrix arrays
  27. var asize = a._size;
  28. // rows and columns
  29. var m = asize[0];
  30. var n = asize[1];
  31. // variables
  32. var i, j, k, J, p, p0, p1;
  33. // workspace size
  34. var s = 4 * n + (ata ? n + m + 1 : 0);
  35. // allocate workspace
  36. var w = []; // (s)
  37. var ancestor = 0; // first n entries
  38. var maxfirst = n; // next n entries
  39. var prevleaf = 2 * n; // next n entries
  40. var first = 3 * n; // next n entries
  41. var head = 4 * n; // next n + 1 entries (used when ata is true)
  42. var next = 5 * n + 1; // last entries in workspace
  43. // clear workspace w[0..s-1]
  44. for (k = 0; k < s; k++) {
  45. w[k] = -1;
  46. }
  47. // allocate result
  48. var colcount = []; // (n)
  49. // AT = A'
  50. var at = transpose(a);
  51. // at arrays
  52. var tindex = at._index;
  53. var tptr = at._ptr;
  54. // find w[first + j]
  55. for (k = 0; k < n; k++) {
  56. j = post[k];
  57. // colcount[j]=1 if j is a leaf
  58. colcount[j] = w[first + j] === -1 ? 1 : 0;
  59. for (; j !== -1 && w[first + j] === -1; j = parent[j]) {
  60. w[first + j] = k;
  61. }
  62. }
  63. // initialize ata if needed
  64. if (ata) {
  65. // invert post
  66. for (k = 0; k < n; k++) {
  67. w[post[k]] = k;
  68. }
  69. // loop rows (columns in AT)
  70. for (i = 0; i < m; i++) {
  71. // values in column i of AT
  72. for (k = n, p0 = tptr[i], p1 = tptr[i + 1], p = p0; p < p1; p++) {
  73. k = Math.min(k, w[tindex[p]]);
  74. }
  75. // place row i in linked list k
  76. w[next + i] = w[head + k];
  77. w[head + k] = i;
  78. }
  79. }
  80. // each node in its own set
  81. for (i = 0; i < n; i++) {
  82. w[ancestor + i] = i;
  83. }
  84. for (k = 0; k < n; k++) {
  85. // j is the kth node in postordered etree
  86. j = post[k];
  87. // check j is not a root
  88. if (parent[j] !== -1) {
  89. colcount[parent[j]]--;
  90. }
  91. // J=j for LL'=A case
  92. for (J = ata ? w[head + k] : j; J !== -1; J = ata ? w[next + J] : -1) {
  93. for (p = tptr[J]; p < tptr[J + 1]; p++) {
  94. i = tindex[p];
  95. var r = csLeaf(i, j, w, first, maxfirst, prevleaf, ancestor);
  96. // check A(i,j) is in skeleton
  97. if (r.jleaf >= 1) {
  98. colcount[j]++;
  99. }
  100. // check account for overlap in q
  101. if (r.jleaf === 2) {
  102. colcount[r.q]--;
  103. }
  104. }
  105. }
  106. if (parent[j] !== -1) {
  107. w[ancestor + j] = parent[j];
  108. }
  109. }
  110. // sum up colcount's of each child
  111. for (j = 0; j < n; j++) {
  112. if (parent[j] !== -1) {
  113. colcount[parent[j]] += colcount[j];
  114. }
  115. }
  116. return colcount;
  117. };
  118. });