1 | import { factory } from '../../../utils/factory.js';
|
2 | import { DimensionError } from '../../../error/DimensionError.js';
|
3 | var name = 'algorithm04';
|
4 | var dependencies = ['typed', 'equalScalar'];
|
5 | export var createAlgorithm04 = /* #__PURE__ */factory(name, dependencies, _ref => {
|
6 | var {
|
7 | typed,
|
8 | equalScalar
|
9 | } = _ref;
|
10 |
|
11 | /**
|
12 | * Iterates over SparseMatrix A and SparseMatrix B nonzero items and invokes the callback function f(Aij, Bij).
|
13 | * Callback function invoked MAX(NNZA, NNZB) times
|
14 | *
|
15 | *
|
16 | * ┌ f(Aij, Bij) ; A(i,j) !== 0 && B(i,j) !== 0
|
17 | * C(i,j) = ┤ A(i,j) ; A(i,j) !== 0
|
18 | * └ B(i,j) ; B(i,j) !== 0
|
19 | *
|
20 | *
|
21 | * @param {Matrix} a The SparseMatrix instance (A)
|
22 | * @param {Matrix} b The SparseMatrix instance (B)
|
23 | * @param {Function} callback The f(Aij,Bij) operation to invoke
|
24 | *
|
25 | * @return {Matrix} SparseMatrix (C)
|
26 | *
|
27 | * see https://github.com/josdejong/mathjs/pull/346#issuecomment-97620294
|
28 | */
|
29 | return function algorithm04(a, b, callback) {
|
30 | // sparse matrix arrays
|
31 | var avalues = a._values;
|
32 | var aindex = a._index;
|
33 | var aptr = a._ptr;
|
34 | var asize = a._size;
|
35 | var adt = a._datatype; // sparse matrix arrays
|
36 |
|
37 | var bvalues = b._values;
|
38 | var bindex = b._index;
|
39 | var bptr = b._ptr;
|
40 | var bsize = b._size;
|
41 | var bdt = b._datatype; // validate dimensions
|
42 |
|
43 | if (asize.length !== bsize.length) {
|
44 | throw new DimensionError(asize.length, bsize.length);
|
45 | } // check rows & columns
|
46 |
|
47 |
|
48 | if (asize[0] !== bsize[0] || asize[1] !== bsize[1]) {
|
49 | throw new RangeError('Dimension mismatch. Matrix A (' + asize + ') must match Matrix B (' + bsize + ')');
|
50 | } // rows & columns
|
51 |
|
52 |
|
53 | var rows = asize[0];
|
54 | var columns = asize[1]; // datatype
|
55 |
|
56 | var dt; // equal signature to use
|
57 |
|
58 | var eq = equalScalar; // zero value
|
59 |
|
60 | var zero = 0; // callback signature to use
|
61 |
|
62 | var cf = callback; // process data types
|
63 |
|
64 | if (typeof adt === 'string' && adt === bdt) {
|
65 | // datatype
|
66 | dt = adt; // find signature that matches (dt, dt)
|
67 |
|
68 | eq = typed.find(equalScalar, [dt, dt]); // convert 0 to the same datatype
|
69 |
|
70 | zero = typed.convert(0, dt); // callback
|
71 |
|
72 | cf = typed.find(callback, [dt, dt]);
|
73 | } // result arrays
|
74 |
|
75 |
|
76 | var cvalues = avalues && bvalues ? [] : undefined;
|
77 | var cindex = [];
|
78 | var cptr = []; // workspace
|
79 |
|
80 | var xa = avalues && bvalues ? [] : undefined;
|
81 | var xb = avalues && bvalues ? [] : undefined; // marks indicating we have a value in x for a given column
|
82 |
|
83 | var wa = [];
|
84 | var wb = []; // vars
|
85 |
|
86 | var i, j, k, k0, k1; // loop columns
|
87 |
|
88 | for (j = 0; j < columns; j++) {
|
89 | // update cptr
|
90 | cptr[j] = cindex.length; // columns mark
|
91 |
|
92 | var mark = j + 1; // loop A(:,j)
|
93 |
|
94 | for (k0 = aptr[j], k1 = aptr[j + 1], k = k0; k < k1; k++) {
|
95 | // row
|
96 | i = aindex[k]; // update c
|
97 |
|
98 | cindex.push(i); // update workspace
|
99 |
|
100 | wa[i] = mark; // check we need to process values
|
101 |
|
102 | if (xa) {
|
103 | xa[i] = avalues[k];
|
104 | }
|
105 | } // loop B(:,j)
|
106 |
|
107 |
|
108 | for (k0 = bptr[j], k1 = bptr[j + 1], k = k0; k < k1; k++) {
|
109 | // row
|
110 | i = bindex[k]; // check row exists in A
|
111 |
|
112 | if (wa[i] === mark) {
|
113 | // update record in xa @ i
|
114 | if (xa) {
|
115 | // invoke callback
|
116 | var v = cf(xa[i], bvalues[k]); // check for zero
|
117 |
|
118 | if (!eq(v, zero)) {
|
119 | // update workspace
|
120 | xa[i] = v;
|
121 | } else {
|
122 | // remove mark (index will be removed later)
|
123 | wa[i] = null;
|
124 | }
|
125 | }
|
126 | } else {
|
127 | // update c
|
128 | cindex.push(i); // update workspace
|
129 |
|
130 | wb[i] = mark; // check we need to process values
|
131 |
|
132 | if (xb) {
|
133 | xb[i] = bvalues[k];
|
134 | }
|
135 | }
|
136 | } // check we need to process values (non pattern matrix)
|
137 |
|
138 |
|
139 | if (xa && xb) {
|
140 | // initialize first index in j
|
141 | k = cptr[j]; // loop index in j
|
142 |
|
143 | while (k < cindex.length) {
|
144 | // row
|
145 | i = cindex[k]; // check workspace has value @ i
|
146 |
|
147 | if (wa[i] === mark) {
|
148 | // push value (Aij != 0 || (Aij != 0 && Bij != 0))
|
149 | cvalues[k] = xa[i]; // increment pointer
|
150 |
|
151 | k++;
|
152 | } else if (wb[i] === mark) {
|
153 | // push value (bij != 0)
|
154 | cvalues[k] = xb[i]; // increment pointer
|
155 |
|
156 | k++;
|
157 | } else {
|
158 | // remove index @ k
|
159 | cindex.splice(k, 1);
|
160 | }
|
161 | }
|
162 | }
|
163 | } // update cptr
|
164 |
|
165 |
|
166 | cptr[columns] = cindex.length; // return sparse matrix
|
167 |
|
168 | return a.createSparseMatrix({
|
169 | values: cvalues,
|
170 | index: cindex,
|
171 | ptr: cptr,
|
172 | size: [rows, columns],
|
173 | datatype: dt
|
174 | });
|
175 | };
|
176 | }); |
\ | No newline at end of file |