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