csSqr.js (5057B)
1 "use strict"; 2 3 Object.defineProperty(exports, "__esModule", { 4 value: true 5 }); 6 exports.createCsSqr = void 0; 7 8 var _csPermute = require("./csPermute.js"); 9 10 var _csPost = require("./csPost.js"); 11 12 var _csEtree = require("./csEtree.js"); 13 14 var _csAmd = require("./csAmd.js"); 15 16 var _csCounts = require("./csCounts.js"); 17 18 var _factory = require("../../../utils/factory.js"); 19 20 var name = 'csSqr'; 21 var dependencies = ['add', 'multiply', 'transpose']; 22 var createCsSqr = /* #__PURE__ */(0, _factory.factory)(name, dependencies, function (_ref) { 23 var add = _ref.add, 24 multiply = _ref.multiply, 25 transpose = _ref.transpose; 26 var csAmd = (0, _csAmd.createCsAmd)({ 27 add: add, 28 multiply: multiply, 29 transpose: transpose 30 }); 31 var csCounts = (0, _csCounts.createCsCounts)({ 32 transpose: transpose 33 }); 34 /** 35 * Symbolic ordering and analysis for QR and LU decompositions. 36 * 37 * @param {Number} order The ordering strategy (see csAmd for more details) 38 * @param {Matrix} a The A matrix 39 * @param {boolean} qr Symbolic ordering and analysis for QR decomposition (true) or 40 * symbolic ordering and analysis for LU decomposition (false) 41 * 42 * @return {Object} The Symbolic ordering and analysis for matrix A 43 * 44 * Reference: http://faculty.cse.tamu.edu/davis/publications.html 45 */ 46 47 return function csSqr(order, a, qr) { 48 // a arrays 49 var aptr = a._ptr; 50 var asize = a._size; // columns 51 52 var n = asize[1]; // vars 53 54 var k; // symbolic analysis result 55 56 var s = {}; // fill-reducing ordering 57 58 s.q = csAmd(order, a); // validate results 59 60 if (order && !s.q) { 61 return null; 62 } // QR symbolic analysis 63 64 65 if (qr) { 66 // apply permutations if needed 67 var c = order ? (0, _csPermute.csPermute)(a, null, s.q, 0) : a; // etree of C'*C, where C=A(:,q) 68 69 s.parent = (0, _csEtree.csEtree)(c, 1); // post order elimination tree 70 71 var post = (0, _csPost.csPost)(s.parent, n); // col counts chol(C'*C) 72 73 s.cp = csCounts(c, s.parent, post, 1); // check we have everything needed to calculate number of nonzero elements 74 75 if (c && s.parent && s.cp && _vcount(c, s)) { 76 // calculate number of nonzero elements 77 for (s.unz = 0, k = 0; k < n; k++) { 78 s.unz += s.cp[k]; 79 } 80 } 81 } else { 82 // for LU factorization only, guess nnz(L) and nnz(U) 83 s.unz = 4 * aptr[n] + n; 84 s.lnz = s.unz; 85 } // return result S 86 87 88 return s; 89 }; 90 /** 91 * Compute nnz(V) = s.lnz, s.pinv, s.leftmost, s.m2 from A and s.parent 92 */ 93 94 function _vcount(a, s) { 95 // a arrays 96 var aptr = a._ptr; 97 var aindex = a._index; 98 var asize = a._size; // rows & columns 99 100 var m = asize[0]; 101 var n = asize[1]; // initialize s arrays 102 103 s.pinv = []; // (m + n) 104 105 s.leftmost = []; // (m) 106 // vars 107 108 var parent = s.parent; 109 var pinv = s.pinv; 110 var leftmost = s.leftmost; // workspace, next: first m entries, head: next n entries, tail: next n entries, nque: next n entries 111 112 var w = []; // (m + 3 * n) 113 114 var next = 0; 115 var head = m; 116 var tail = m + n; 117 var nque = m + 2 * n; // vars 118 119 var i, k, p, p0, p1; // initialize w 120 121 for (k = 0; k < n; k++) { 122 // queue k is empty 123 w[head + k] = -1; 124 w[tail + k] = -1; 125 w[nque + k] = 0; 126 } // initialize row arrays 127 128 129 for (i = 0; i < m; i++) { 130 leftmost[i] = -1; 131 } // loop columns backwards 132 133 134 for (k = n - 1; k >= 0; k--) { 135 // values & index for column k 136 for (p0 = aptr[k], p1 = aptr[k + 1], p = p0; p < p1; p++) { 137 // leftmost[i] = min(find(A(i,:))) 138 leftmost[aindex[p]] = k; 139 } 140 } // scan rows in reverse order 141 142 143 for (i = m - 1; i >= 0; i--) { 144 // row i is not yet ordered 145 pinv[i] = -1; 146 k = leftmost[i]; // check row i is empty 147 148 if (k === -1) { 149 continue; 150 } // first row in queue k 151 152 153 if (w[nque + k]++ === 0) { 154 w[tail + k] = i; 155 } // put i at head of queue k 156 157 158 w[next + i] = w[head + k]; 159 w[head + k] = i; 160 } 161 162 s.lnz = 0; 163 s.m2 = m; // find row permutation and nnz(V) 164 165 for (k = 0; k < n; k++) { 166 // remove row i from queue k 167 i = w[head + k]; // count V(k,k) as nonzero 168 169 s.lnz++; // add a fictitious row 170 171 if (i < 0) { 172 i = s.m2++; 173 } // associate row i with V(:,k) 174 175 176 pinv[i] = k; // skip if V(k+1:m,k) is empty 177 178 if (--nque[k] <= 0) { 179 continue; 180 } // nque[k] is nnz (V(k+1:m,k)) 181 182 183 s.lnz += w[nque + k]; // move all rows to parent of k 184 185 var pa = parent[k]; 186 187 if (pa !== -1) { 188 if (w[nque + pa] === 0) { 189 w[tail + pa] = w[tail + k]; 190 } 191 192 w[next + w[tail + k]] = w[head + pa]; 193 w[head + pa] = w[next + i]; 194 w[nque + pa] += w[nque + k]; 195 } 196 } 197 198 for (i = 0; i < m; i++) { 199 if (pinv[i] < 0) { 200 pinv[i] = k++; 201 } 202 } 203 204 return true; 205 } 206 }); 207 exports.createCsSqr = createCsSqr;