|
5164
|
1 /* |
|
|
2 |
|
|
3 Copyright (C) 2004 David Bateman |
|
7016
|
4 Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 Andy Adler |
|
|
5 |
|
|
6 This file is part of Octave. |
|
5164
|
7 |
|
|
8 Octave is free software; you can redistribute it and/or modify it |
|
|
9 under the terms of the GNU General Public License as published by the |
|
7016
|
10 Free Software Foundation; either version 3 of the License, or (at your |
|
|
11 option) any later version. |
|
5164
|
12 |
|
|
13 Octave is distributed in the hope that it will be useful, but WITHOUT |
|
|
14 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or |
|
|
15 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License |
|
|
16 for more details. |
|
|
17 |
|
|
18 You should have received a copy of the GNU General Public License |
|
7016
|
19 along with Octave; see the file COPYING. If not, see |
|
|
20 <http://www.gnu.org/licenses/>. |
|
5164
|
21 |
|
|
22 */ |
|
|
23 |
|
|
24 #if !defined (octave_sparse_op_defs_h) |
|
|
25 #define octave_sparse_op_defs_h 1 |
|
|
26 |
|
|
27 #include "Array-util.h" |
|
6221
|
28 #include "mx-ops.h" |
|
5164
|
29 |
|
6708
|
30 #define SPARSE_BIN_OP_DECL(R, OP, X, Y, API) \ |
|
|
31 extern API R OP (const X&, const Y&) |
|
5164
|
32 |
|
6708
|
33 #define SPARSE_CMP_OP_DECL(OP, X, Y, API) \ |
|
|
34 extern API SparseBoolMatrix OP (const X&, const Y&) |
|
5164
|
35 |
|
6708
|
36 #define SPARSE_BOOL_OP_DECL(OP, X, Y, API) \ |
|
|
37 extern API SparseBoolMatrix OP (const X&, const Y&) |
|
5164
|
38 |
|
|
39 // matrix by scalar operations. |
|
|
40 |
|
6708
|
41 #define SPARSE_SMS_BIN_OP_DECLS(R1, R2, M, S, API) \ |
|
|
42 SPARSE_BIN_OP_DECL (R1, operator +, M, S, API); \ |
|
|
43 SPARSE_BIN_OP_DECL (R1, operator -, M, S, API); \ |
|
|
44 SPARSE_BIN_OP_DECL (R2, operator *, M, S, API); \ |
|
|
45 SPARSE_BIN_OP_DECL (R2, operator /, M, S, API); |
|
5164
|
46 |
|
|
47 #define SPARSE_SMS_BIN_OP_1(R, F, OP, M, S) \ |
|
|
48 R \ |
|
|
49 F (const M& m, const S& s) \ |
|
|
50 { \ |
|
5275
|
51 octave_idx_type nr = m.rows (); \ |
|
|
52 octave_idx_type nc = m.cols (); \ |
|
5164
|
53 \ |
|
|
54 R r (nr, nc, (0.0 OP s)); \ |
|
|
55 \ |
|
5275
|
56 for (octave_idx_type j = 0; j < nc; j++) \ |
|
|
57 for (octave_idx_type i = m.cidx (j); i < m.cidx (j+1); i++) \ |
|
5164
|
58 r.elem (m.ridx (i), j) = m.data (i) OP s; \ |
|
|
59 return r; \ |
|
|
60 } |
|
|
61 |
|
|
62 #define SPARSE_SMS_BIN_OP_2(R, F, OP, M, S) \ |
|
|
63 R \ |
|
|
64 F (const M& m, const S& s) \ |
|
|
65 { \ |
|
5275
|
66 octave_idx_type nr = m.rows (); \ |
|
|
67 octave_idx_type nc = m.cols (); \ |
|
5681
|
68 octave_idx_type nz = m.nnz (); \ |
|
5164
|
69 \ |
|
|
70 R r (nr, nc, nz); \ |
|
|
71 \ |
|
5275
|
72 for (octave_idx_type i = 0; i < nz; i++) \ |
|
5164
|
73 { \ |
|
|
74 r.data(i) = m.data(i) OP s; \ |
|
|
75 r.ridx(i) = m.ridx(i); \ |
|
|
76 } \ |
|
5275
|
77 for (octave_idx_type i = 0; i < nc + 1; i++) \ |
|
5164
|
78 r.cidx(i) = m.cidx(i); \ |
|
|
79 \ |
|
|
80 r.maybe_compress (true); \ |
|
|
81 return r; \ |
|
|
82 } |
|
|
83 |
|
|
84 #define SPARSE_SMS_BIN_OPS(R1, R2, M, S) \ |
|
|
85 SPARSE_SMS_BIN_OP_1 (R1, operator +, +, M, S) \ |
|
|
86 SPARSE_SMS_BIN_OP_1 (R1, operator -, -, M, S) \ |
|
|
87 SPARSE_SMS_BIN_OP_2 (R2, operator *, *, M, S) \ |
|
|
88 SPARSE_SMS_BIN_OP_2 (R2, operator /, /, M, S) |
|
|
89 |
|
6708
|
90 #define SPARSE_SMS_CMP_OP_DECLS(M, S, API) \ |
|
|
91 SPARSE_CMP_OP_DECL (mx_el_lt, M, S, API); \ |
|
|
92 SPARSE_CMP_OP_DECL (mx_el_le, M, S, API); \ |
|
|
93 SPARSE_CMP_OP_DECL (mx_el_ge, M, S, API); \ |
|
|
94 SPARSE_CMP_OP_DECL (mx_el_gt, M, S, API); \ |
|
|
95 SPARSE_CMP_OP_DECL (mx_el_eq, M, S, API); \ |
|
|
96 SPARSE_CMP_OP_DECL (mx_el_ne, M, S, API); |
|
5164
|
97 |
|
6708
|
98 #define SPARSE_SMS_EQNE_OP_DECLS(M, S, API) \ |
|
|
99 SPARSE_CMP_OP_DECL (mx_el_eq, M, S, API); \ |
|
|
100 SPARSE_CMP_OP_DECL (mx_el_ne, M, S, API); |
|
5164
|
101 |
|
|
102 #define SPARSE_SMS_CMP_OP(F, OP, M, MZ, MC, S, SZ, SC) \ |
|
|
103 SparseBoolMatrix \ |
|
|
104 F (const M& m, const S& s) \ |
|
|
105 { \ |
|
|
106 /* Count num of non-zero elements */ \ |
|
5275
|
107 octave_idx_type nel = 0; \ |
|
5681
|
108 octave_idx_type nz = m.nnz (); \ |
|
5164
|
109 if (MC (MZ) OP SC (s)) \ |
|
|
110 nel += m.numel() - nz; \ |
|
5275
|
111 for (octave_idx_type i = 0; i < nz; i++) \ |
|
5164
|
112 if (MC (m.data (i)) OP SC (s)) \ |
|
|
113 nel++; \ |
|
|
114 \ |
|
5275
|
115 octave_idx_type nr = m.rows (); \ |
|
|
116 octave_idx_type nc = m.cols (); \ |
|
5164
|
117 SparseBoolMatrix r (nr, nc, nel); \ |
|
|
118 \ |
|
|
119 if (nr > 0 && nc > 0) \ |
|
|
120 { \ |
|
|
121 if (MC (MZ) OP SC (s)) \ |
|
|
122 { \ |
|
5275
|
123 octave_idx_type ii = 0; \ |
|
5164
|
124 r.cidx (0) = 0; \ |
|
5275
|
125 for (octave_idx_type j = 0; j < nc; j++) \ |
|
5164
|
126 { \ |
|
5275
|
127 for (octave_idx_type i = 0; i < nr; i++) \ |
|
5164
|
128 { \ |
|
|
129 bool el = MC (m.elem(i, j)) OP SC (s); \ |
|
|
130 if (el) \ |
|
|
131 { \ |
|
|
132 r.data(ii) = el; \ |
|
|
133 r.ridx(ii++) = i; \ |
|
|
134 } \ |
|
|
135 } \ |
|
|
136 r.cidx(j+1) = ii; \ |
|
|
137 } \ |
|
|
138 } \ |
|
|
139 else \ |
|
|
140 { \ |
|
5275
|
141 octave_idx_type ii = 0; \ |
|
5164
|
142 r.cidx (0) = 0; \ |
|
5275
|
143 for (octave_idx_type j = 0; j < nc; j++) \ |
|
5164
|
144 { \ |
|
5275
|
145 for (octave_idx_type i = m.cidx(j); i < m.cidx(j+1); i++) \ |
|
5164
|
146 { \ |
|
|
147 bool el = MC (m.data(i)) OP SC (s); \ |
|
|
148 if (el) \ |
|
|
149 { \ |
|
|
150 r.data(ii) = el; \ |
|
|
151 r.ridx(ii++) = m.ridx(i); \ |
|
|
152 } \ |
|
|
153 } \ |
|
|
154 r.cidx(j+1) = ii; \ |
|
|
155 } \ |
|
|
156 } \ |
|
|
157 } \ |
|
|
158 return r; \ |
|
|
159 } |
|
|
160 |
|
|
161 #define SPARSE_SMS_CMP_OPS(M, MZ, CM, S, SZ, CS) \ |
|
|
162 SPARSE_SMS_CMP_OP (mx_el_lt, <, M, MZ, CM, S, SZ, CS) \ |
|
|
163 SPARSE_SMS_CMP_OP (mx_el_le, <=, M, MZ, CM, S, SZ, CS) \ |
|
|
164 SPARSE_SMS_CMP_OP (mx_el_ge, >=, M, MZ, CM, S, SZ, CS) \ |
|
|
165 SPARSE_SMS_CMP_OP (mx_el_gt, >, M, MZ, CM, S, SZ, CS) \ |
|
|
166 SPARSE_SMS_CMP_OP (mx_el_eq, ==, M, MZ, , S, SZ, ) \ |
|
|
167 SPARSE_SMS_CMP_OP (mx_el_ne, !=, M, MZ, , S, SZ, ) |
|
|
168 |
|
|
169 #define SPARSE_SMS_EQNE_OPS(M, MZ, CM, S, SZ, CS) \ |
|
|
170 SPARSE_SMS_CMP_OP (mx_el_eq, ==, M, MZ, , S, SZ, ) \ |
|
|
171 SPARSE_SMS_CMP_OP (mx_el_ne, !=, M, MZ, , S, SZ, ) |
|
|
172 |
|
6708
|
173 #define SPARSE_SMS_BOOL_OP_DECLS(M, S, API) \ |
|
|
174 SPARSE_BOOL_OP_DECL (mx_el_and, M, S, API); \ |
|
|
175 SPARSE_BOOL_OP_DECL (mx_el_or, M, S, API); |
|
5164
|
176 |
|
|
177 #define SPARSE_SMS_BOOL_OP(F, OP, M, S, LHS_ZERO, RHS_ZERO) \ |
|
|
178 SparseBoolMatrix \ |
|
|
179 F (const M& m, const S& s) \ |
|
|
180 { \ |
|
|
181 /* Count num of non-zero elements */ \ |
|
5275
|
182 octave_idx_type nel = 0; \ |
|
5681
|
183 octave_idx_type nz = m.nnz (); \ |
|
5164
|
184 if (LHS_ZERO OP (s != RHS_ZERO)) \ |
|
|
185 nel += m.numel() - nz; \ |
|
5275
|
186 for (octave_idx_type i = 0; i < nz; i++) \ |
|
5164
|
187 if ((m.data(i) != LHS_ZERO) OP (s != RHS_ZERO))\ |
|
|
188 nel++; \ |
|
|
189 \ |
|
5275
|
190 octave_idx_type nr = m.rows (); \ |
|
|
191 octave_idx_type nc = m.cols (); \ |
|
5164
|
192 SparseBoolMatrix r (nr, nc, nel); \ |
|
|
193 \ |
|
|
194 if (nr > 0 && nc > 0) \ |
|
|
195 { \ |
|
|
196 if (LHS_ZERO OP (s != RHS_ZERO)) \ |
|
|
197 { \ |
|
5275
|
198 octave_idx_type ii = 0; \ |
|
5164
|
199 r.cidx (0) = 0; \ |
|
5275
|
200 for (octave_idx_type j = 0; j < nc; j++) \ |
|
5164
|
201 { \ |
|
5275
|
202 for (octave_idx_type i = 0; i < nr; i++) \ |
|
5164
|
203 { \ |
|
|
204 bool el = (m.elem(i, j) != LHS_ZERO) OP (s != RHS_ZERO); \ |
|
|
205 if (el) \ |
|
|
206 { \ |
|
|
207 r.data(ii) = el; \ |
|
|
208 r.ridx(ii++) = i; \ |
|
|
209 } \ |
|
|
210 } \ |
|
|
211 r.cidx(j+1) = ii; \ |
|
|
212 } \ |
|
|
213 } \ |
|
|
214 else \ |
|
|
215 { \ |
|
5275
|
216 octave_idx_type ii = 0; \ |
|
5164
|
217 r.cidx (0) = 0; \ |
|
5275
|
218 for (octave_idx_type j = 0; j < nc; j++) \ |
|
5164
|
219 { \ |
|
5275
|
220 for (octave_idx_type i = m.cidx(j); i < m.cidx(j+1); i++) \ |
|
5164
|
221 { \ |
|
|
222 bool el = (m.data(i) != LHS_ZERO) OP (s != RHS_ZERO); \ |
|
|
223 if (el) \ |
|
|
224 { \ |
|
|
225 r.data(ii) = el; \ |
|
|
226 r.ridx(ii++) = m.ridx(i); \ |
|
|
227 } \ |
|
|
228 } \ |
|
|
229 r.cidx(j+1) = ii; \ |
|
|
230 } \ |
|
|
231 } \ |
|
|
232 } \ |
|
|
233 return r; \ |
|
|
234 } |
|
|
235 |
|
|
236 #define SPARSE_SMS_BOOL_OPS2(M, S, LHS_ZERO, RHS_ZERO) \ |
|
|
237 SPARSE_SMS_BOOL_OP (mx_el_and, &&, M, S, LHS_ZERO, RHS_ZERO) \ |
|
|
238 SPARSE_SMS_BOOL_OP (mx_el_or, ||, M, S, LHS_ZERO, RHS_ZERO) |
|
|
239 |
|
|
240 #define SPARSE_SMS_BOOL_OPS(M, S, ZERO) \ |
|
|
241 SPARSE_SMS_BOOL_OPS2(M, S, ZERO, ZERO) |
|
|
242 |
|
6708
|
243 #define SPARSE_SMS_OP_DECLS(R1, R2, M, S, API) \ |
|
|
244 SPARSE_SMS_BIN_OP_DECLS (R1, R2, M, S, API) \ |
|
|
245 SPARSE_SMS_CMP_OP_DECLS (M, S, API) \ |
|
|
246 SPARSE_SMS_BOOL_OP_DECLS (M, S, API) |
|
5164
|
247 |
|
|
248 // scalar by matrix operations. |
|
|
249 |
|
6708
|
250 #define SPARSE_SSM_BIN_OP_DECLS(R1, R2, S, M, API) \ |
|
|
251 SPARSE_BIN_OP_DECL (R1, operator +, S, M, API); \ |
|
|
252 SPARSE_BIN_OP_DECL (R1, operator -, S, M, API); \ |
|
|
253 SPARSE_BIN_OP_DECL (R2, operator *, S, M, API); \ |
|
|
254 SPARSE_BIN_OP_DECL (R2, operator /, S, M, API); |
|
5164
|
255 |
|
|
256 #define SPARSE_SSM_BIN_OP_1(R, F, OP, S, M) \ |
|
|
257 R \ |
|
|
258 F (const S& s, const M& m) \ |
|
|
259 { \ |
|
5275
|
260 octave_idx_type nr = m.rows (); \ |
|
|
261 octave_idx_type nc = m.cols (); \ |
|
5164
|
262 \ |
|
|
263 R r (nr, nc, (s OP 0.0)); \ |
|
|
264 \ |
|
5275
|
265 for (octave_idx_type j = 0; j < nc; j++) \ |
|
|
266 for (octave_idx_type i = m.cidx (j); i < m.cidx (j+1); i++) \ |
|
5164
|
267 r.elem (m.ridx (i), j) = s OP m.data (i); \ |
|
|
268 \ |
|
|
269 return r; \ |
|
|
270 } |
|
|
271 |
|
|
272 #define SPARSE_SSM_BIN_OP_2(R, F, OP, S, M) \ |
|
|
273 R \ |
|
|
274 F (const S& s, const M& m) \ |
|
|
275 { \ |
|
5275
|
276 octave_idx_type nr = m.rows (); \ |
|
|
277 octave_idx_type nc = m.cols (); \ |
|
5681
|
278 octave_idx_type nz = m.nnz (); \ |
|
5164
|
279 \ |
|
|
280 R r (nr, nc, nz); \ |
|
|
281 \ |
|
5275
|
282 for (octave_idx_type i = 0; i < nz; i++) \ |
|
5164
|
283 { \ |
|
|
284 r.data(i) = s OP m.data(i); \ |
|
|
285 r.ridx(i) = m.ridx(i); \ |
|
|
286 } \ |
|
5275
|
287 for (octave_idx_type i = 0; i < nc + 1; i++) \ |
|
5164
|
288 r.cidx(i) = m.cidx(i); \ |
|
|
289 \ |
|
|
290 r.maybe_compress(true); \ |
|
|
291 return r; \ |
|
|
292 } |
|
|
293 |
|
|
294 #define SPARSE_SSM_BIN_OPS(R1, R2, S, M) \ |
|
|
295 SPARSE_SSM_BIN_OP_1 (R1, operator +, +, S, M) \ |
|
|
296 SPARSE_SSM_BIN_OP_1 (R1, operator -, -, S, M) \ |
|
|
297 SPARSE_SSM_BIN_OP_2 (R2, operator *, *, S, M) \ |
|
|
298 SPARSE_SSM_BIN_OP_2 (R2, operator /, /, S, M) |
|
|
299 |
|
6708
|
300 #define SPARSE_SSM_CMP_OP_DECLS(S, M, API) \ |
|
|
301 SPARSE_CMP_OP_DECL (mx_el_lt, S, M, API); \ |
|
|
302 SPARSE_CMP_OP_DECL (mx_el_le, S, M, API); \ |
|
|
303 SPARSE_CMP_OP_DECL (mx_el_ge, S, M, API); \ |
|
|
304 SPARSE_CMP_OP_DECL (mx_el_gt, S, M, API); \ |
|
|
305 SPARSE_CMP_OP_DECL (mx_el_eq, S, M, API); \ |
|
|
306 SPARSE_CMP_OP_DECL (mx_el_ne, S, M, API); |
|
5164
|
307 |
|
6708
|
308 #define SPARSE_SSM_EQNE_OP_DECLS(S, M, API) \ |
|
|
309 SPARSE_CMP_OP_DECL (mx_el_eq, S, M, API); \ |
|
|
310 SPARSE_CMP_OP_DECL (mx_el_ne, S, M, API); |
|
5164
|
311 |
|
|
312 #define SPARSE_SSM_CMP_OP(F, OP, S, SZ, SC, M, MZ, MC) \ |
|
|
313 SparseBoolMatrix \ |
|
|
314 F (const S& s, const M& m) \ |
|
|
315 { \ |
|
|
316 /* Count num of non-zero elements */ \ |
|
5275
|
317 octave_idx_type nel = 0; \ |
|
5681
|
318 octave_idx_type nz = m.nnz (); \ |
|
5164
|
319 if (SC (s) OP MC (MZ)) \ |
|
|
320 nel += m.numel() - nz; \ |
|
5275
|
321 for (octave_idx_type i = 0; i < nz; i++) \ |
|
5164
|
322 if (SC (s) OP MC (m.data (i))) \ |
|
|
323 nel++; \ |
|
|
324 \ |
|
5275
|
325 octave_idx_type nr = m.rows (); \ |
|
|
326 octave_idx_type nc = m.cols (); \ |
|
5164
|
327 SparseBoolMatrix r (nr, nc, nel); \ |
|
|
328 \ |
|
|
329 if (nr > 0 && nc > 0) \ |
|
|
330 { \ |
|
|
331 if (SC (s) OP MC (MZ))\ |
|
|
332 { \ |
|
5275
|
333 octave_idx_type ii = 0; \ |
|
5164
|
334 r.cidx (0) = 0; \ |
|
5275
|
335 for (octave_idx_type j = 0; j < nc; j++) \ |
|
5164
|
336 { \ |
|
5275
|
337 for (octave_idx_type i = 0; i < nr; i++) \ |
|
5164
|
338 { \ |
|
|
339 bool el = SC (s) OP MC (m.elem(i, j)); \ |
|
|
340 if (el) \ |
|
|
341 { \ |
|
|
342 r.data(ii) = el; \ |
|
|
343 r.ridx(ii++) = i; \ |
|
|
344 } \ |
|
|
345 } \ |
|
|
346 r.cidx(j+1) = ii; \ |
|
|
347 } \ |
|
|
348 } \ |
|
|
349 else \ |
|
|
350 { \ |
|
5275
|
351 octave_idx_type ii = 0; \ |
|
5164
|
352 r.cidx (0) = 0; \ |
|
5275
|
353 for (octave_idx_type j = 0; j < nc; j++) \ |
|
5164
|
354 { \ |
|
5275
|
355 for (octave_idx_type i = m.cidx(j); i < m.cidx(j+1); i++) \ |
|
5164
|
356 { \ |
|
|
357 bool el = SC (s) OP MC (m.data(i)); \ |
|
|
358 if (el) \ |
|
|
359 { \ |
|
|
360 r.data(ii) = el; \ |
|
|
361 r.ridx(ii++) = m.ridx(i); \ |
|
|
362 } \ |
|
|
363 } \ |
|
|
364 r.cidx(j+1) = ii; \ |
|
|
365 } \ |
|
|
366 } \ |
|
|
367 } \ |
|
|
368 return r; \ |
|
|
369 } |
|
|
370 |
|
|
371 #define SPARSE_SSM_CMP_OPS(S, SZ, SC, M, MZ, MC) \ |
|
|
372 SPARSE_SSM_CMP_OP (mx_el_lt, <, S, SZ, SC, M, MZ, MC) \ |
|
|
373 SPARSE_SSM_CMP_OP (mx_el_le, <=, S, SZ, SC, M, MZ, MC) \ |
|
|
374 SPARSE_SSM_CMP_OP (mx_el_ge, >=, S, SZ, SC, M, MZ, MC) \ |
|
|
375 SPARSE_SSM_CMP_OP (mx_el_gt, >, S, SZ, SC, M, MZ, MC) \ |
|
|
376 SPARSE_SSM_CMP_OP (mx_el_eq, ==, S, SZ, , M, MZ, ) \ |
|
|
377 SPARSE_SSM_CMP_OP (mx_el_ne, !=, S, SZ, , M, MZ, ) |
|
|
378 |
|
|
379 #define SPARSE_SSM_EQNE_OPS(S, SZ, SC, M, MZ, MC) \ |
|
|
380 SPARSE_SSM_CMP_OP (mx_el_eq, ==, S, SZ, , M, MZ, ) \ |
|
|
381 SPARSE_SSM_CMP_OP (mx_el_ne, !=, S, SZ, , M, MZ, ) |
|
|
382 |
|
6708
|
383 #define SPARSE_SSM_BOOL_OP_DECLS(S, M, API) \ |
|
|
384 SPARSE_BOOL_OP_DECL (mx_el_and, S, M, API); \ |
|
|
385 SPARSE_BOOL_OP_DECL (mx_el_or, S, M, API); \ |
|
5164
|
386 |
|
|
387 #define SPARSE_SSM_BOOL_OP(F, OP, S, M, LHS_ZERO, RHS_ZERO) \ |
|
|
388 SparseBoolMatrix \ |
|
|
389 F (const S& s, const M& m) \ |
|
|
390 { \ |
|
|
391 /* Count num of non-zero elements */ \ |
|
5275
|
392 octave_idx_type nel = 0; \ |
|
5681
|
393 octave_idx_type nz = m.nnz (); \ |
|
5164
|
394 if ((s != LHS_ZERO) OP RHS_ZERO) \ |
|
|
395 nel += m.numel() - nz; \ |
|
5275
|
396 for (octave_idx_type i = 0; i < nz; i++) \ |
|
5164
|
397 if ((s != LHS_ZERO) OP m.data(i) != RHS_ZERO) \ |
|
|
398 nel++; \ |
|
|
399 \ |
|
5275
|
400 octave_idx_type nr = m.rows (); \ |
|
|
401 octave_idx_type nc = m.cols (); \ |
|
5164
|
402 SparseBoolMatrix r (nr, nc, nel); \ |
|
|
403 \ |
|
|
404 if (nr > 0 && nc > 0) \ |
|
|
405 { \ |
|
|
406 if ((s != LHS_ZERO) OP RHS_ZERO) \ |
|
|
407 { \ |
|
5275
|
408 octave_idx_type ii = 0; \ |
|
5164
|
409 r.cidx (0) = 0; \ |
|
5275
|
410 for (octave_idx_type j = 0; j < nc; j++) \ |
|
5164
|
411 { \ |
|
5275
|
412 for (octave_idx_type i = 0; i < nr; i++) \ |
|
5164
|
413 { \ |
|
|
414 bool el = (s != LHS_ZERO) OP (m.elem(i, j) != RHS_ZERO); \ |
|
|
415 if (el) \ |
|
|
416 { \ |
|
|
417 r.data(ii) = el; \ |
|
|
418 r.ridx(ii++) = i; \ |
|
|
419 } \ |
|
|
420 } \ |
|
|
421 r.cidx(j+1) = ii; \ |
|
|
422 } \ |
|
|
423 } \ |
|
|
424 else \ |
|
|
425 { \ |
|
5275
|
426 octave_idx_type ii = 0; \ |
|
5164
|
427 r.cidx (0) = 0; \ |
|
5275
|
428 for (octave_idx_type j = 0; j < nc; j++) \ |
|
5164
|
429 { \ |
|
5275
|
430 for (octave_idx_type i = m.cidx(j); i < m.cidx(j+1); i++) \ |
|
5164
|
431 { \ |
|
|
432 bool el = (s != LHS_ZERO) OP (m.data(i) != RHS_ZERO); \ |
|
|
433 if (el) \ |
|
|
434 { \ |
|
|
435 r.data(ii) = el; \ |
|
|
436 r.ridx(ii++) = m.ridx(i); \ |
|
|
437 } \ |
|
|
438 } \ |
|
|
439 r.cidx(j+1) = ii; \ |
|
|
440 } \ |
|
|
441 } \ |
|
|
442 } \ |
|
|
443 return r; \ |
|
|
444 } |
|
|
445 |
|
|
446 #define SPARSE_SSM_BOOL_OPS2(S, M, LHS_ZERO, RHS_ZERO) \ |
|
|
447 SPARSE_SSM_BOOL_OP (mx_el_and, &&, S, M, LHS_ZERO, RHS_ZERO) \ |
|
|
448 SPARSE_SSM_BOOL_OP (mx_el_or, ||, S, M, LHS_ZERO, RHS_ZERO) |
|
|
449 |
|
|
450 #define SPARSE_SSM_BOOL_OPS(S, M, ZERO) \ |
|
|
451 SPARSE_SSM_BOOL_OPS2(S, M, ZERO, ZERO) |
|
|
452 |
|
6708
|
453 #define SPARSE_SSM_OP_DECLS(R1, R2, S, M, API) \ |
|
|
454 SPARSE_SSM_BIN_OP_DECLS (R1, R2, S, M, API) \ |
|
|
455 SPARSE_SSM_CMP_OP_DECLS (S, M, API) \ |
|
|
456 SPARSE_SSM_BOOL_OP_DECLS (S, M, API) \ |
|
5164
|
457 |
|
|
458 // matrix by matrix operations. |
|
|
459 |
|
6708
|
460 #define SPARSE_SMSM_BIN_OP_DECLS(R1, R2, M1, M2, API) \ |
|
|
461 SPARSE_BIN_OP_DECL (R1, operator +, M1, M2, API); \ |
|
|
462 SPARSE_BIN_OP_DECL (R1, operator -, M1, M2, API); \ |
|
|
463 SPARSE_BIN_OP_DECL (R2, product, M1, M2, API); \ |
|
|
464 SPARSE_BIN_OP_DECL (R2, quotient, M1, M2, API); |
|
5164
|
465 |
|
|
466 #define SPARSE_SMSM_BIN_OP_1(R, F, OP, M1, M2) \ |
|
|
467 R \ |
|
|
468 F (const M1& m1, const M2& m2) \ |
|
|
469 { \ |
|
|
470 R r; \ |
|
|
471 \ |
|
5275
|
472 octave_idx_type m1_nr = m1.rows (); \ |
|
|
473 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
474 \ |
|
5275
|
475 octave_idx_type m2_nr = m2.rows (); \ |
|
|
476 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
477 \ |
|
6221
|
478 if (m1_nr == 1 && m1_nc == 1) \ |
|
|
479 { \ |
|
|
480 if (m1.elem(0,0) == 0.) \ |
|
|
481 r = R (m2); \ |
|
|
482 else \ |
|
|
483 { \ |
|
|
484 r = R (m2_nr, m2_nc, m1.data(0) OP 0.); \ |
|
|
485 \ |
|
|
486 for (octave_idx_type j = 0 ; j < m2_nc ; j++) \ |
|
|
487 { \ |
|
|
488 OCTAVE_QUIT; \ |
|
|
489 octave_idx_type idxj = j * m2_nr; \ |
|
|
490 for (octave_idx_type i = m2.cidx(j) ; i < m2.cidx(j+1) ; i++) \ |
|
|
491 { \ |
|
|
492 OCTAVE_QUIT; \ |
|
|
493 r.data(idxj + m2.ridx(i)) = m1.data(0) OP m2.data(i); \ |
|
|
494 } \ |
|
|
495 } \ |
|
|
496 r.maybe_compress (); \ |
|
|
497 } \ |
|
|
498 } \ |
|
|
499 else if (m2_nr == 1 && m2_nc == 1) \ |
|
|
500 { \ |
|
|
501 if (m2.elem(0,0) == 0.) \ |
|
|
502 r = R (m1); \ |
|
|
503 else \ |
|
|
504 { \ |
|
|
505 r = R (m1_nr, m1_nc, 0. OP m2.data(0)); \ |
|
|
506 \ |
|
|
507 for (octave_idx_type j = 0 ; j < m1_nc ; j++) \ |
|
|
508 { \ |
|
|
509 OCTAVE_QUIT; \ |
|
|
510 octave_idx_type idxj = j * m1_nr; \ |
|
|
511 for (octave_idx_type i = m1.cidx(j) ; i < m1.cidx(j+1) ; i++) \ |
|
|
512 { \ |
|
|
513 OCTAVE_QUIT; \ |
|
|
514 r.data(idxj + m1.ridx(i)) = m1.data(i) OP m2.data(0); \ |
|
|
515 } \ |
|
|
516 } \ |
|
|
517 r.maybe_compress (); \ |
|
|
518 } \ |
|
|
519 } \ |
|
|
520 else if (m1_nr != m2_nr || m1_nc != m2_nc) \ |
|
5164
|
521 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
522 else \ |
|
|
523 { \ |
|
5681
|
524 r = R (m1_nr, m1_nc, (m1.nnz () + m2.nnz ())); \ |
|
5164
|
525 \ |
|
5275
|
526 octave_idx_type jx = 0; \ |
|
5164
|
527 r.cidx (0) = 0; \ |
|
5275
|
528 for (octave_idx_type i = 0 ; i < m1_nc ; i++) \ |
|
5164
|
529 { \ |
|
5275
|
530 octave_idx_type ja = m1.cidx(i); \ |
|
|
531 octave_idx_type ja_max = m1.cidx(i+1); \ |
|
5164
|
532 bool ja_lt_max= ja < ja_max; \ |
|
|
533 \ |
|
5275
|
534 octave_idx_type jb = m2.cidx(i); \ |
|
|
535 octave_idx_type jb_max = m2.cidx(i+1); \ |
|
5164
|
536 bool jb_lt_max = jb < jb_max; \ |
|
|
537 \ |
|
|
538 while (ja_lt_max || jb_lt_max ) \ |
|
|
539 { \ |
|
|
540 OCTAVE_QUIT; \ |
|
|
541 if ((! jb_lt_max) || \ |
|
|
542 (ja_lt_max && (m1.ridx(ja) < m2.ridx(jb)))) \ |
|
|
543 { \ |
|
|
544 r.ridx(jx) = m1.ridx(ja); \ |
|
|
545 r.data(jx) = m1.data(ja) OP 0.; \ |
|
|
546 jx++; \ |
|
|
547 ja++; \ |
|
|
548 ja_lt_max= ja < ja_max; \ |
|
|
549 } \ |
|
|
550 else if (( !ja_lt_max ) || \ |
|
|
551 (jb_lt_max && (m2.ridx(jb) < m1.ridx(ja)) ) ) \ |
|
|
552 { \ |
|
|
553 r.ridx(jx) = m2.ridx(jb); \ |
|
|
554 r.data(jx) = 0. OP m2.data(jb); \ |
|
|
555 jx++; \ |
|
|
556 jb++; \ |
|
|
557 jb_lt_max= jb < jb_max; \ |
|
|
558 } \ |
|
|
559 else \ |
|
|
560 { \ |
|
|
561 if ((m1.data(ja) OP m2.data(jb)) != 0.) \ |
|
|
562 { \ |
|
|
563 r.data(jx) = m1.data(ja) OP m2.data(jb); \ |
|
|
564 r.ridx(jx) = m1.ridx(ja); \ |
|
|
565 jx++; \ |
|
|
566 } \ |
|
|
567 ja++; \ |
|
|
568 ja_lt_max= ja < ja_max; \ |
|
|
569 jb++; \ |
|
|
570 jb_lt_max= jb < jb_max; \ |
|
|
571 } \ |
|
|
572 } \ |
|
|
573 r.cidx(i+1) = jx; \ |
|
|
574 } \ |
|
|
575 \ |
|
|
576 r.maybe_compress (); \ |
|
|
577 } \ |
|
|
578 \ |
|
|
579 return r; \ |
|
|
580 } |
|
|
581 |
|
|
582 #define SPARSE_SMSM_BIN_OP_2(R, F, OP, M1, M2) \ |
|
|
583 R \ |
|
|
584 F (const M1& m1, const M2& m2) \ |
|
|
585 { \ |
|
|
586 R r; \ |
|
|
587 \ |
|
5275
|
588 octave_idx_type m1_nr = m1.rows (); \ |
|
|
589 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
590 \ |
|
5275
|
591 octave_idx_type m2_nr = m2.rows (); \ |
|
|
592 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
593 \ |
|
6221
|
594 if (m1_nr == 1 && m1_nc == 1) \ |
|
|
595 { \ |
|
|
596 if (m1.elem(0,0) == 0.) \ |
|
|
597 r = R (m2_nr, m2_nc); \ |
|
|
598 else \ |
|
|
599 { \ |
|
|
600 r = R (m2); \ |
|
|
601 octave_idx_type m2_nnz = m2.nnz(); \ |
|
|
602 \ |
|
|
603 for (octave_idx_type i = 0 ; i < m2_nnz ; i++) \ |
|
|
604 { \ |
|
|
605 OCTAVE_QUIT; \ |
|
|
606 r.data (i) = m1.data(0) OP r.data(i); \ |
|
|
607 } \ |
|
|
608 r.maybe_compress (); \ |
|
|
609 } \ |
|
|
610 } \ |
|
|
611 else if (m2_nr == 1 && m2_nc == 1) \ |
|
|
612 { \ |
|
|
613 if (m2.elem(0,0) == 0.) \ |
|
|
614 r = R (m1_nr, m1_nc); \ |
|
|
615 else \ |
|
|
616 { \ |
|
|
617 r = R (m1); \ |
|
|
618 octave_idx_type m1_nnz = m1.nnz(); \ |
|
|
619 \ |
|
|
620 for (octave_idx_type i = 0 ; i < m1_nnz ; i++) \ |
|
|
621 { \ |
|
|
622 OCTAVE_QUIT; \ |
|
|
623 r.data (i) = r.data(i) OP m2.data(0); \ |
|
|
624 } \ |
|
|
625 r.maybe_compress (); \ |
|
|
626 } \ |
|
|
627 } \ |
|
|
628 else if (m1_nr != m2_nr || m1_nc != m2_nc) \ |
|
5164
|
629 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
630 else \ |
|
|
631 { \ |
|
5681
|
632 r = R (m1_nr, m1_nc, (m1.nnz () > m2.nnz () ? m1.nnz () : m2.nnz ())); \ |
|
5164
|
633 \ |
|
5275
|
634 octave_idx_type jx = 0; \ |
|
5164
|
635 r.cidx (0) = 0; \ |
|
5275
|
636 for (octave_idx_type i = 0 ; i < m1_nc ; i++) \ |
|
5164
|
637 { \ |
|
5275
|
638 octave_idx_type ja = m1.cidx(i); \ |
|
|
639 octave_idx_type ja_max = m1.cidx(i+1); \ |
|
5164
|
640 bool ja_lt_max= ja < ja_max; \ |
|
|
641 \ |
|
5275
|
642 octave_idx_type jb = m2.cidx(i); \ |
|
|
643 octave_idx_type jb_max = m2.cidx(i+1); \ |
|
5164
|
644 bool jb_lt_max = jb < jb_max; \ |
|
|
645 \ |
|
|
646 while (ja_lt_max || jb_lt_max ) \ |
|
|
647 { \ |
|
|
648 OCTAVE_QUIT; \ |
|
|
649 if ((! jb_lt_max) || \ |
|
|
650 (ja_lt_max && (m1.ridx(ja) < m2.ridx(jb)))) \ |
|
|
651 { \ |
|
|
652 ja++; ja_lt_max= ja < ja_max; \ |
|
|
653 } \ |
|
|
654 else if (( !ja_lt_max ) || \ |
|
|
655 (jb_lt_max && (m2.ridx(jb) < m1.ridx(ja)) ) ) \ |
|
|
656 { \ |
|
|
657 jb++; jb_lt_max= jb < jb_max; \ |
|
|
658 } \ |
|
|
659 else \ |
|
|
660 { \ |
|
|
661 if ((m1.data(ja) OP m2.data(jb)) != 0.) \ |
|
|
662 { \ |
|
|
663 r.data(jx) = m1.data(ja) OP m2.data(jb); \ |
|
|
664 r.ridx(jx) = m1.ridx(ja); \ |
|
|
665 jx++; \ |
|
|
666 } \ |
|
|
667 ja++; ja_lt_max= ja < ja_max; \ |
|
|
668 jb++; jb_lt_max= jb < jb_max; \ |
|
|
669 } \ |
|
|
670 } \ |
|
|
671 r.cidx(i+1) = jx; \ |
|
|
672 } \ |
|
|
673 \ |
|
|
674 r.maybe_compress (); \ |
|
|
675 } \ |
|
|
676 \ |
|
|
677 return r; \ |
|
|
678 } |
|
|
679 |
|
|
680 #define SPARSE_SMSM_BIN_OP_3(R, F, OP, M1, M2) \ |
|
|
681 R \ |
|
|
682 F (const M1& m1, const M2& m2) \ |
|
|
683 { \ |
|
|
684 R r; \ |
|
|
685 \ |
|
5275
|
686 octave_idx_type m1_nr = m1.rows (); \ |
|
|
687 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
688 \ |
|
5275
|
689 octave_idx_type m2_nr = m2.rows (); \ |
|
|
690 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
691 \ |
|
6221
|
692 if (m1_nr == 1 && m1_nc == 1) \ |
|
|
693 { \ |
|
|
694 if ((m1.elem (0,0) OP Complex()) == Complex()) \ |
|
|
695 { \ |
|
|
696 octave_idx_type m2_nnz = m2.nnz(); \ |
|
|
697 r = R (m2); \ |
|
|
698 for (octave_idx_type i = 0 ; i < m2_nnz ; i++) \ |
|
|
699 r.data (i) = m1.elem(0,0) OP r.data(i); \ |
|
|
700 r.maybe_compress (); \ |
|
|
701 } \ |
|
|
702 else \ |
|
|
703 { \ |
|
|
704 r = R (m2_nr, m2_nc, m1.elem(0,0) OP Complex ()); \ |
|
|
705 for (octave_idx_type j = 0 ; j < m2_nc ; j++) \ |
|
|
706 { \ |
|
|
707 OCTAVE_QUIT; \ |
|
|
708 octave_idx_type idxj = j * m2_nr; \ |
|
|
709 for (octave_idx_type i = m2.cidx(j) ; i < m2.cidx(j+1) ; i++) \ |
|
|
710 { \ |
|
|
711 OCTAVE_QUIT; \ |
|
|
712 r.data(idxj + m2.ridx(i)) = m1.elem(0,0) OP m2.data(i); \ |
|
|
713 } \ |
|
|
714 } \ |
|
|
715 r.maybe_compress (); \ |
|
|
716 } \ |
|
|
717 } \ |
|
|
718 else if (m2_nr == 1 && m2_nc == 1) \ |
|
|
719 { \ |
|
|
720 if ((Complex() OP m1.elem (0,0)) == Complex()) \ |
|
|
721 { \ |
|
|
722 octave_idx_type m1_nnz = m1.nnz(); \ |
|
|
723 r = R (m1); \ |
|
|
724 for (octave_idx_type i = 0 ; i < m1_nnz ; i++) \ |
|
|
725 r.data (i) = r.data(i) OP m2.elem(0,0); \ |
|
|
726 r.maybe_compress (); \ |
|
|
727 } \ |
|
|
728 else \ |
|
|
729 { \ |
|
|
730 r = R (m1_nr, m1_nc, Complex() OP m2.elem(0,0)); \ |
|
|
731 for (octave_idx_type j = 0 ; j < m1_nc ; j++) \ |
|
|
732 { \ |
|
|
733 OCTAVE_QUIT; \ |
|
|
734 octave_idx_type idxj = j * m1_nr; \ |
|
|
735 for (octave_idx_type i = m1.cidx(j) ; i < m1.cidx(j+1) ; i++) \ |
|
|
736 { \ |
|
|
737 OCTAVE_QUIT; \ |
|
|
738 r.data(idxj + m1.ridx(i)) = m1.data(i) OP m2.elem(0,0); \ |
|
|
739 } \ |
|
|
740 } \ |
|
|
741 r.maybe_compress (); \ |
|
|
742 } \ |
|
|
743 } \ |
|
|
744 else if (m1_nr != m2_nr || m1_nc != m2_nc) \ |
|
5164
|
745 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
746 else \ |
|
|
747 { \ |
|
|
748 \ |
|
5775
|
749 /* FIXME Kludge... Always double/Complex, so Complex () */ \ |
|
5164
|
750 r = R (m1_nr, m1_nc, (Complex () OP Complex ())); \ |
|
|
751 \ |
|
5275
|
752 for (octave_idx_type i = 0 ; i < m1_nc ; i++) \ |
|
5164
|
753 { \ |
|
5275
|
754 octave_idx_type ja = m1.cidx(i); \ |
|
|
755 octave_idx_type ja_max = m1.cidx(i+1); \ |
|
5164
|
756 bool ja_lt_max= ja < ja_max; \ |
|
|
757 \ |
|
5275
|
758 octave_idx_type jb = m2.cidx(i); \ |
|
|
759 octave_idx_type jb_max = m2.cidx(i+1); \ |
|
5164
|
760 bool jb_lt_max = jb < jb_max; \ |
|
|
761 \ |
|
|
762 while (ja_lt_max || jb_lt_max ) \ |
|
|
763 { \ |
|
|
764 OCTAVE_QUIT; \ |
|
|
765 if ((! jb_lt_max) || \ |
|
|
766 (ja_lt_max && (m1.ridx(ja) < m2.ridx(jb)))) \ |
|
|
767 { \ |
|
|
768 /* keep those kludges coming */ \ |
|
|
769 r.elem(m1.ridx(ja),i) = m1.data(ja) OP Complex (); \ |
|
|
770 ja++; \ |
|
|
771 ja_lt_max= ja < ja_max; \ |
|
|
772 } \ |
|
|
773 else if (( !ja_lt_max ) || \ |
|
|
774 (jb_lt_max && (m2.ridx(jb) < m1.ridx(ja)) ) ) \ |
|
|
775 { \ |
|
|
776 /* keep those kludges coming */ \ |
|
|
777 r.elem(m2.ridx(jb),i) = Complex () OP m2.data(jb); \ |
|
|
778 jb++; \ |
|
|
779 jb_lt_max= jb < jb_max; \ |
|
|
780 } \ |
|
|
781 else \ |
|
|
782 { \ |
|
|
783 r.elem(m1.ridx(ja),i) = m1.data(ja) OP m2.data(jb); \ |
|
|
784 ja++; \ |
|
|
785 ja_lt_max= ja < ja_max; \ |
|
|
786 jb++; \ |
|
|
787 jb_lt_max= jb < jb_max; \ |
|
|
788 } \ |
|
|
789 } \ |
|
|
790 } \ |
|
|
791 r.maybe_compress (true); \ |
|
|
792 } \ |
|
|
793 \ |
|
|
794 return r; \ |
|
|
795 } |
|
|
796 |
|
|
797 // Note that SM ./ SM needs to take into account the NaN and Inf values |
|
|
798 // implied by the division by zero. |
|
5775
|
799 // FIXME Are the NaNs double(NaN) or Complex(NaN,Nan) in the complex |
|
5164
|
800 // case? |
|
|
801 #define SPARSE_SMSM_BIN_OPS(R1, R2, M1, M2) \ |
|
|
802 SPARSE_SMSM_BIN_OP_1 (R1, operator +, +, M1, M2) \ |
|
|
803 SPARSE_SMSM_BIN_OP_1 (R1, operator -, -, M1, M2) \ |
|
|
804 SPARSE_SMSM_BIN_OP_2 (R2, product, *, M1, M2) \ |
|
|
805 SPARSE_SMSM_BIN_OP_3 (R2, quotient, /, M1, M2) |
|
|
806 |
|
6708
|
807 #define SPARSE_SMSM_CMP_OP_DECLS(M1, M2, API) \ |
|
|
808 SPARSE_CMP_OP_DECL (mx_el_lt, M1, M2, API); \ |
|
|
809 SPARSE_CMP_OP_DECL (mx_el_le, M1, M2, API); \ |
|
|
810 SPARSE_CMP_OP_DECL (mx_el_ge, M1, M2, API); \ |
|
|
811 SPARSE_CMP_OP_DECL (mx_el_gt, M1, M2, API); \ |
|
|
812 SPARSE_CMP_OP_DECL (mx_el_eq, M1, M2, API); \ |
|
|
813 SPARSE_CMP_OP_DECL (mx_el_ne, M1, M2, API); |
|
5164
|
814 |
|
6708
|
815 #define SPARSE_SMSM_EQNE_OP_DECLS(M1, M2, API) \ |
|
|
816 SPARSE_CMP_OP_DECL (mx_el_eq, M1, M2, API); \ |
|
|
817 SPARSE_CMP_OP_DECL (mx_el_ne, M1, M2, API); |
|
5164
|
818 |
|
|
819 #define SPARSE_SMSM_CMP_OP(F, OP, M1, C1, M2, C2) \ |
|
|
820 SparseBoolMatrix \ |
|
|
821 F (const M1& m1, const M2& m2) \ |
|
|
822 { \ |
|
|
823 SparseBoolMatrix r; \ |
|
|
824 \ |
|
5275
|
825 octave_idx_type m1_nr = m1.rows (); \ |
|
|
826 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
827 \ |
|
5275
|
828 octave_idx_type m2_nr = m2.rows (); \ |
|
|
829 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
830 \ |
|
6221
|
831 if (m1_nr == 1 && m1_nc == 1) \ |
|
|
832 { \ |
|
|
833 extern OCTAVE_API SparseBoolMatrix F (const double&, const M2&); \ |
|
|
834 extern OCTAVE_API SparseBoolMatrix F (const Complex&, const M2&); \ |
|
|
835 r = F (m1.elem(0,0), m2); \ |
|
|
836 } \ |
|
|
837 else if (m2_nr == 1 && m2_nc == 1) \ |
|
|
838 { \ |
|
|
839 extern OCTAVE_API SparseBoolMatrix F (const M1&, const double&); \ |
|
|
840 extern OCTAVE_API SparseBoolMatrix F (const M1&, const Complex&); \ |
|
|
841 r = F (m1, m2.elem(0,0)); \ |
|
|
842 } \ |
|
|
843 else if (m1_nr == m2_nr && m1_nc == m2_nc) \ |
|
5164
|
844 { \ |
|
|
845 if (m1_nr != 0 || m1_nc != 0) \ |
|
|
846 { \ |
|
|
847 /* Count num of non-zero elements */ \ |
|
5275
|
848 octave_idx_type nel = 0; \ |
|
|
849 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
|
850 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
851 if (C1 (m1.elem(i, j)) OP C2 (m2.elem(i, j))) \ |
|
|
852 nel++; \ |
|
|
853 \ |
|
|
854 r = SparseBoolMatrix (m1_nr, m1_nc, nel); \ |
|
|
855 \ |
|
5275
|
856 octave_idx_type ii = 0; \ |
|
5164
|
857 r.cidx (0) = 0; \ |
|
5275
|
858 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
5164
|
859 { \ |
|
5275
|
860 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
861 { \ |
|
|
862 bool el = C1 (m1.elem(i, j)) OP C2 (m2.elem(i, j)); \ |
|
|
863 if (el) \ |
|
|
864 { \ |
|
|
865 r.data(ii) = el; \ |
|
|
866 r.ridx(ii++) = i; \ |
|
|
867 } \ |
|
|
868 } \ |
|
|
869 r.cidx(j+1) = ii; \ |
|
|
870 } \ |
|
|
871 } \ |
|
|
872 } \ |
|
|
873 else \ |
|
|
874 { \ |
|
|
875 if ((m1_nr != 0 || m1_nc != 0) && (m2_nr != 0 || m2_nc != 0)) \ |
|
|
876 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
877 } \ |
|
|
878 return r; \ |
|
|
879 } |
|
|
880 |
|
|
881 #define SPARSE_SMSM_CMP_OPS(M1, Z1, C1, M2, Z2, C2) \ |
|
|
882 SPARSE_SMSM_CMP_OP (mx_el_lt, <, M1, C1, M2, C2) \ |
|
|
883 SPARSE_SMSM_CMP_OP (mx_el_le, <=, M1, C1, M2, C2) \ |
|
|
884 SPARSE_SMSM_CMP_OP (mx_el_ge, >=, M1, C1, M2, C2) \ |
|
|
885 SPARSE_SMSM_CMP_OP (mx_el_gt, >, M1, C1, M2, C2) \ |
|
|
886 SPARSE_SMSM_CMP_OP (mx_el_eq, ==, M1, , M2, ) \ |
|
|
887 SPARSE_SMSM_CMP_OP (mx_el_ne, !=, M1, , M2, ) |
|
|
888 |
|
|
889 #define SPARSE_SMSM_EQNE_OPS(M1, Z1, C1, M2, Z2, C2) \ |
|
|
890 SPARSE_SMSM_CMP_OP (mx_el_eq, ==, M1, , M2, ) \ |
|
|
891 SPARSE_SMSM_CMP_OP (mx_el_ne, !=, M1, , M2, ) |
|
|
892 |
|
6708
|
893 #define SPARSE_SMSM_BOOL_OP_DECLS(M1, M2, API) \ |
|
|
894 SPARSE_BOOL_OP_DECL (mx_el_and, M1, M2, API); \ |
|
|
895 SPARSE_BOOL_OP_DECL (mx_el_or, M1, M2, API); |
|
5164
|
896 |
|
|
897 #define SPARSE_SMSM_BOOL_OP(F, OP, M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
898 SparseBoolMatrix \ |
|
|
899 F (const M1& m1, const M2& m2) \ |
|
|
900 { \ |
|
|
901 SparseBoolMatrix r; \ |
|
|
902 \ |
|
5275
|
903 octave_idx_type m1_nr = m1.rows (); \ |
|
|
904 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
905 \ |
|
5275
|
906 octave_idx_type m2_nr = m2.rows (); \ |
|
|
907 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
908 \ |
|
6221
|
909 if (m1_nr == 1 && m1_nc == 1) \ |
|
|
910 { \ |
|
|
911 extern OCTAVE_API SparseBoolMatrix F (const double&, const M2&); \ |
|
|
912 extern OCTAVE_API SparseBoolMatrix F (const Complex&, const M2&); \ |
|
|
913 r = F (m1.elem(0,0), m2); \ |
|
|
914 } \ |
|
|
915 else if (m2_nr == 1 && m2_nc == 1) \ |
|
|
916 { \ |
|
|
917 extern OCTAVE_API SparseBoolMatrix F (const M1&, const double&); \ |
|
|
918 extern OCTAVE_API SparseBoolMatrix F (const M1&, const Complex&); \ |
|
|
919 r = F (m1, m2.elem(0,0)); \ |
|
|
920 } \ |
|
|
921 else if (m1_nr == m2_nr && m1_nc == m2_nc) \ |
|
5164
|
922 { \ |
|
|
923 if (m1_nr != 0 || m1_nc != 0) \ |
|
|
924 { \ |
|
|
925 /* Count num of non-zero elements */ \ |
|
5275
|
926 octave_idx_type nel = 0; \ |
|
|
927 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
|
928 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
929 if ((m1.elem(i, j) != LHS_ZERO) \ |
|
|
930 OP (m2.elem(i, j) != RHS_ZERO)) \ |
|
|
931 nel++; \ |
|
|
932 \ |
|
|
933 r = SparseBoolMatrix (m1_nr, m1_nc, nel); \ |
|
|
934 \ |
|
5275
|
935 octave_idx_type ii = 0; \ |
|
5164
|
936 r.cidx (0) = 0; \ |
|
5275
|
937 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
5164
|
938 { \ |
|
5275
|
939 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
940 { \ |
|
|
941 bool el = (m1.elem(i, j) != LHS_ZERO) \ |
|
|
942 OP (m2.elem(i, j) != RHS_ZERO); \ |
|
|
943 if (el) \ |
|
|
944 { \ |
|
|
945 r.data(ii) = el; \ |
|
|
946 r.ridx(ii++) = i; \ |
|
|
947 } \ |
|
|
948 } \ |
|
|
949 r.cidx(j+1) = ii; \ |
|
|
950 } \ |
|
|
951 } \ |
|
|
952 } \ |
|
|
953 else \ |
|
|
954 { \ |
|
|
955 if ((m1_nr != 0 || m1_nc != 0) && (m2_nr != 0 || m2_nc != 0)) \ |
|
|
956 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
957 } \ |
|
|
958 return r; \ |
|
|
959 } |
|
|
960 |
|
|
961 #define SPARSE_SMSM_BOOL_OPS2(M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
962 SPARSE_SMSM_BOOL_OP (mx_el_and, &&, M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
963 SPARSE_SMSM_BOOL_OP (mx_el_or, ||, M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
964 |
|
|
965 #define SPARSE_SMSM_BOOL_OPS(M1, M2, ZERO) \ |
|
|
966 SPARSE_SMSM_BOOL_OPS2(M1, M2, ZERO, ZERO) |
|
|
967 |
|
6708
|
968 #define SPARSE_SMSM_OP_DECLS(R1, R2, M1, M2, API) \ |
|
|
969 SPARSE_SMSM_BIN_OP_DECLS (R1, R2, M1, M2, API) \ |
|
|
970 SPARSE_SMSM_CMP_OP_DECLS (M1, M2, API) \ |
|
|
971 SPARSE_SMSM_BOOL_OP_DECLS (M1, M2, API) |
|
5164
|
972 |
|
|
973 // matrix by matrix operations. |
|
|
974 |
|
6708
|
975 #define SPARSE_MSM_BIN_OP_DECLS(R1, R2, M1, M2, API) \ |
|
|
976 SPARSE_BIN_OP_DECL (R1, operator +, M1, M2, API); \ |
|
|
977 SPARSE_BIN_OP_DECL (R1, operator -, M1, M2, API); \ |
|
|
978 SPARSE_BIN_OP_DECL (R2, product, M1, M2, API); \ |
|
|
979 SPARSE_BIN_OP_DECL (R2, quotient, M1, M2, API); |
|
5164
|
980 |
|
|
981 #define SPARSE_MSM_BIN_OP_1(R, F, OP, M1, M2) \ |
|
|
982 R \ |
|
|
983 F (const M1& m1, const M2& m2) \ |
|
|
984 { \ |
|
|
985 R r; \ |
|
|
986 \ |
|
5275
|
987 octave_idx_type m1_nr = m1.rows (); \ |
|
|
988 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
989 \ |
|
5275
|
990 octave_idx_type m2_nr = m2.rows (); \ |
|
|
991 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
992 \ |
|
6221
|
993 if (m2_nr == 1 && m2_nc == 1) \ |
|
|
994 r = R (m1 OP m2.elem(0,0)); \ |
|
|
995 else if (m1_nr != m2_nr || m1_nc != m2_nc) \ |
|
5164
|
996 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
997 else \ |
|
|
998 { \ |
|
|
999 r = R (m1_nr, m1_nc); \ |
|
|
1000 \ |
|
5275
|
1001 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
|
1002 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1003 r.elem (i, j) = m1.elem (i, j) OP m2.elem (i, j); \ |
|
|
1004 } \ |
|
|
1005 return r; \ |
|
|
1006 } |
|
|
1007 |
|
|
1008 #define SPARSE_MSM_BIN_OP_2(R, F, OP, M1, M2, ZERO) \ |
|
|
1009 R \ |
|
|
1010 F (const M1& m1, const M2& m2) \ |
|
|
1011 { \ |
|
|
1012 R r; \ |
|
|
1013 \ |
|
5275
|
1014 octave_idx_type m1_nr = m1.rows (); \ |
|
|
1015 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
1016 \ |
|
5275
|
1017 octave_idx_type m2_nr = m2.rows (); \ |
|
|
1018 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
1019 \ |
|
6221
|
1020 if (m2_nr == 1 && m2_nc == 1) \ |
|
|
1021 r = R (m1 OP m2.elem(0,0)); \ |
|
|
1022 else if (m1_nr != m2_nr || m1_nc != m2_nc) \ |
|
5164
|
1023 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
1024 else \ |
|
|
1025 { \ |
|
|
1026 /* Count num of non-zero elements */ \ |
|
5275
|
1027 octave_idx_type nel = 0; \ |
|
|
1028 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
|
1029 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1030 if ((m1.elem(i, j) OP m2.elem(i, j)) != ZERO) \ |
|
|
1031 nel++; \ |
|
|
1032 \ |
|
|
1033 r = R (m1_nr, m1_nc, nel); \ |
|
|
1034 \ |
|
5275
|
1035 octave_idx_type ii = 0; \ |
|
5164
|
1036 r.cidx (0) = 0; \ |
|
5275
|
1037 for (octave_idx_type j = 0 ; j < m1_nc ; j++) \ |
|
5164
|
1038 { \ |
|
5275
|
1039 for (octave_idx_type i = 0 ; i < m1_nr ; i++) \ |
|
5164
|
1040 { \ |
|
|
1041 if ((m1.elem(i, j) OP m2.elem(i, j)) != ZERO) \ |
|
|
1042 { \ |
|
|
1043 r.data (ii) = m1.elem(i, j) OP m2.elem(i,j); \ |
|
|
1044 r.ridx (ii++) = i; \ |
|
|
1045 } \ |
|
|
1046 } \ |
|
|
1047 r.cidx(j+1) = ii; \ |
|
|
1048 } \ |
|
|
1049 } \ |
|
|
1050 \ |
|
|
1051 return r; \ |
|
|
1052 } |
|
|
1053 |
|
5775
|
1054 // FIXME Pass a specific ZERO value |
|
5164
|
1055 #define SPARSE_MSM_BIN_OPS(R1, R2, M1, M2) \ |
|
|
1056 SPARSE_MSM_BIN_OP_1 (R1, operator +, +, M1, M2) \ |
|
|
1057 SPARSE_MSM_BIN_OP_1 (R1, operator -, -, M1, M2) \ |
|
|
1058 SPARSE_MSM_BIN_OP_2 (R2, product, *, M1, M2, 0.0) \ |
|
|
1059 SPARSE_MSM_BIN_OP_2 (R2, quotient, /, M1, M2, 0.0) |
|
|
1060 |
|
6708
|
1061 #define SPARSE_MSM_CMP_OP_DECLS(M1, M2, API) \ |
|
|
1062 SPARSE_CMP_OP_DECL (mx_el_lt, M1, M2, API); \ |
|
|
1063 SPARSE_CMP_OP_DECL (mx_el_le, M1, M2, API); \ |
|
|
1064 SPARSE_CMP_OP_DECL (mx_el_ge, M1, M2, API); \ |
|
|
1065 SPARSE_CMP_OP_DECL (mx_el_gt, M1, M2, API); \ |
|
|
1066 SPARSE_CMP_OP_DECL (mx_el_eq, M1, M2, API); \ |
|
|
1067 SPARSE_CMP_OP_DECL (mx_el_ne, M1, M2, API); |
|
5164
|
1068 |
|
6708
|
1069 #define SPARSE_MSM_EQNE_OP_DECLS(M1, M2, API) \ |
|
|
1070 SPARSE_CMP_OP_DECL (mx_el_eq, M1, M2, API); \ |
|
|
1071 SPARSE_CMP_OP_DECL (mx_el_ne, M1, M2, API); |
|
5164
|
1072 |
|
|
1073 #define SPARSE_MSM_CMP_OP(F, OP, M1, C1, M2, C2) \ |
|
|
1074 SparseBoolMatrix \ |
|
|
1075 F (const M1& m1, const M2& m2) \ |
|
|
1076 { \ |
|
|
1077 SparseBoolMatrix r; \ |
|
|
1078 \ |
|
5275
|
1079 octave_idx_type m1_nr = m1.rows (); \ |
|
|
1080 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
1081 \ |
|
5275
|
1082 octave_idx_type m2_nr = m2.rows (); \ |
|
|
1083 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
1084 \ |
|
6221
|
1085 if (m2_nr == 1 && m2_nc == 1) \ |
|
|
1086 r = SparseBoolMatrix (F (m1, m2.elem(0,0))); \ |
|
|
1087 else if (m1_nr == m2_nr && m1_nc == m2_nc) \ |
|
5164
|
1088 { \ |
|
|
1089 if (m1_nr != 0 || m1_nc != 0) \ |
|
|
1090 { \ |
|
|
1091 /* Count num of non-zero elements */ \ |
|
5275
|
1092 octave_idx_type nel = 0; \ |
|
|
1093 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
|
1094 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1095 if (C1 (m1.elem(i, j)) OP C2 (m2.elem(i, j))) \ |
|
|
1096 nel++; \ |
|
|
1097 \ |
|
|
1098 r = SparseBoolMatrix (m1_nr, m1_nc, nel); \ |
|
|
1099 \ |
|
5275
|
1100 octave_idx_type ii = 0; \ |
|
5164
|
1101 r.cidx (0) = 0; \ |
|
5275
|
1102 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
5164
|
1103 { \ |
|
5275
|
1104 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1105 { \ |
|
|
1106 bool el = C1 (m1.elem(i, j)) OP C2 (m2.elem(i, j)); \ |
|
|
1107 if (el) \ |
|
|
1108 { \ |
|
|
1109 r.data(ii) = el; \ |
|
|
1110 r.ridx(ii++) = i; \ |
|
|
1111 } \ |
|
|
1112 } \ |
|
|
1113 r.cidx(j+1) = ii; \ |
|
|
1114 } \ |
|
|
1115 } \ |
|
|
1116 } \ |
|
|
1117 else \ |
|
|
1118 { \ |
|
|
1119 if ((m1_nr != 0 || m1_nc != 0) && (m2_nr != 0 || m2_nc != 0)) \ |
|
|
1120 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
1121 } \ |
|
|
1122 return r; \ |
|
|
1123 } |
|
|
1124 |
|
|
1125 #define SPARSE_MSM_CMP_OPS(M1, Z1, C1, M2, Z2, C2) \ |
|
|
1126 SPARSE_MSM_CMP_OP (mx_el_lt, <, M1, C1, M2, C2) \ |
|
|
1127 SPARSE_MSM_CMP_OP (mx_el_le, <=, M1, C1, M2, C2) \ |
|
|
1128 SPARSE_MSM_CMP_OP (mx_el_ge, >=, M1, C1, M2, C2) \ |
|
|
1129 SPARSE_MSM_CMP_OP (mx_el_gt, >, M1, C1, M2, C2) \ |
|
|
1130 SPARSE_MSM_CMP_OP (mx_el_eq, ==, M1, , M2, ) \ |
|
|
1131 SPARSE_MSM_CMP_OP (mx_el_ne, !=, M1, , M2, ) |
|
|
1132 |
|
|
1133 #define SPARSE_MSM_EQNE_OPS(M1, Z1, C1, M2, Z2, C2) \ |
|
|
1134 SPARSE_MSM_CMP_OP (mx_el_eq, ==, M1, , M2, ) \ |
|
|
1135 SPARSE_MSM_CMP_OP (mx_el_ne, !=, M1, , M2, ) |
|
|
1136 |
|
6708
|
1137 #define SPARSE_MSM_BOOL_OP_DECLS(M1, M2, API) \ |
|
|
1138 SPARSE_BOOL_OP_DECL (mx_el_and, M1, M2, API); \ |
|
|
1139 SPARSE_BOOL_OP_DECL (mx_el_or, M1, M2, API); |
|
5164
|
1140 |
|
|
1141 #define SPARSE_MSM_BOOL_OP(F, OP, M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
1142 SparseBoolMatrix \ |
|
|
1143 F (const M1& m1, const M2& m2) \ |
|
|
1144 { \ |
|
|
1145 SparseBoolMatrix r; \ |
|
|
1146 \ |
|
5275
|
1147 octave_idx_type m1_nr = m1.rows (); \ |
|
|
1148 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
1149 \ |
|
5275
|
1150 octave_idx_type m2_nr = m2.rows (); \ |
|
|
1151 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
1152 \ |
|
6221
|
1153 if (m2_nr == 1 && m2_nc == 1) \ |
|
|
1154 r = SparseBoolMatrix (F (m1, m2.elem(0,0))); \ |
|
|
1155 else if (m1_nr == m2_nr && m1_nc == m2_nc) \ |
|
5164
|
1156 { \ |
|
|
1157 if (m1_nr != 0 || m1_nc != 0) \ |
|
|
1158 { \ |
|
|
1159 /* Count num of non-zero elements */ \ |
|
5275
|
1160 octave_idx_type nel = 0; \ |
|
|
1161 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
|
1162 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1163 if ((m1.elem(i, j) != LHS_ZERO) \ |
|
|
1164 OP (m2.elem(i, j) != RHS_ZERO)) \ |
|
|
1165 nel++; \ |
|
|
1166 \ |
|
|
1167 r = SparseBoolMatrix (m1_nr, m1_nc, nel); \ |
|
|
1168 \ |
|
5275
|
1169 octave_idx_type ii = 0; \ |
|
5164
|
1170 r.cidx (0) = 0; \ |
|
5275
|
1171 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
5164
|
1172 { \ |
|
5275
|
1173 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1174 { \ |
|
|
1175 bool el = (m1.elem(i, j) != LHS_ZERO) \ |
|
|
1176 OP (m2.elem(i, j) != RHS_ZERO); \ |
|
|
1177 if (el) \ |
|
|
1178 { \ |
|
|
1179 r.data(ii) = el; \ |
|
|
1180 r.ridx(ii++) = i; \ |
|
|
1181 } \ |
|
|
1182 } \ |
|
|
1183 r.cidx(j+1) = ii; \ |
|
|
1184 } \ |
|
|
1185 } \ |
|
|
1186 } \ |
|
|
1187 else \ |
|
|
1188 { \ |
|
|
1189 if ((m1_nr != 0 || m1_nc != 0) && (m2_nr != 0 || m2_nc != 0)) \ |
|
|
1190 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
1191 } \ |
|
|
1192 return r; \ |
|
|
1193 } |
|
|
1194 |
|
|
1195 #define SPARSE_MSM_BOOL_OPS2(M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
1196 SPARSE_MSM_BOOL_OP (mx_el_and, &&, M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
1197 SPARSE_MSM_BOOL_OP (mx_el_or, ||, M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
1198 |
|
|
1199 #define SPARSE_MSM_BOOL_OPS(M1, M2, ZERO) \ |
|
|
1200 SPARSE_MSM_BOOL_OPS2(M1, M2, ZERO, ZERO) |
|
|
1201 |
|
6708
|
1202 #define SPARSE_MSM_OP_DECLS(R1, R2, M1, M2, API) \ |
|
|
1203 SPARSE_MSM_BIN_OP_DECLS (R1, R2, M1, M2, API) \ |
|
|
1204 SPARSE_MSM_CMP_OP_DECLS (M1, M2, API) \ |
|
|
1205 SPARSE_MSM_BOOL_OP_DECLS (M1, M2, API) |
|
5164
|
1206 |
|
|
1207 // matrix by matrix operations. |
|
|
1208 |
|
6708
|
1209 #define SPARSE_SMM_BIN_OP_DECLS(R1, R2, M1, M2, API) \ |
|
|
1210 SPARSE_BIN_OP_DECL (R1, operator +, M1, M2, API); \ |
|
|
1211 SPARSE_BIN_OP_DECL (R1, operator -, M1, M2, API); \ |
|
|
1212 SPARSE_BIN_OP_DECL (R2, product, M1, M2, API); \ |
|
|
1213 SPARSE_BIN_OP_DECL (R2, quotient, M1, M2, API); |
|
5164
|
1214 |
|
|
1215 #define SPARSE_SMM_BIN_OP_1(R, F, OP, M1, M2) \ |
|
|
1216 R \ |
|
|
1217 F (const M1& m1, const M2& m2) \ |
|
|
1218 { \ |
|
|
1219 R r; \ |
|
|
1220 \ |
|
5275
|
1221 octave_idx_type m1_nr = m1.rows (); \ |
|
|
1222 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
1223 \ |
|
5275
|
1224 octave_idx_type m2_nr = m2.rows (); \ |
|
|
1225 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
1226 \ |
|
6221
|
1227 if (m1_nr == 1 && m1_nc == 1) \ |
|
|
1228 r = R (m1.elem(0,0) OP m2); \ |
|
|
1229 else if (m1_nr != m2_nr || m1_nc != m2_nc) \ |
|
5164
|
1230 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
1231 else \ |
|
|
1232 { \ |
|
|
1233 r = R (m1_nr, m1_nc); \ |
|
|
1234 \ |
|
5275
|
1235 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
|
1236 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1237 r.elem (i, j) = m1.elem (i, j) OP m2.elem (i, j); \ |
|
|
1238 } \ |
|
|
1239 return r; \ |
|
|
1240 } |
|
|
1241 |
|
|
1242 #define SPARSE_SMM_BIN_OP_2(R, F, OP, M1, M2, ZERO) \ |
|
|
1243 R \ |
|
|
1244 F (const M1& m1, const M2& m2) \ |
|
|
1245 { \ |
|
|
1246 R r; \ |
|
|
1247 \ |
|
5275
|
1248 octave_idx_type m1_nr = m1.rows (); \ |
|
|
1249 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
1250 \ |
|
5275
|
1251 octave_idx_type m2_nr = m2.rows (); \ |
|
|
1252 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
1253 \ |
|
6221
|
1254 if (m1_nr == 1 && m1_nc == 1) \ |
|
|
1255 r = R (m1.elem(0,0) OP m2); \ |
|
|
1256 else if (m1_nr != m2_nr || m1_nc != m2_nc) \ |
|
5164
|
1257 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
1258 else \ |
|
|
1259 { \ |
|
|
1260 /* Count num of non-zero elements */ \ |
|
5275
|
1261 octave_idx_type nel = 0; \ |
|
|
1262 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
|
1263 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1264 if ((m1.elem(i, j) OP m2.elem(i, j)) != ZERO) \ |
|
|
1265 nel++; \ |
|
|
1266 \ |
|
|
1267 r = R (m1_nr, m1_nc, nel); \ |
|
|
1268 \ |
|
5275
|
1269 octave_idx_type ii = 0; \ |
|
5164
|
1270 r.cidx (0) = 0; \ |
|
5275
|
1271 for (octave_idx_type j = 0 ; j < m1_nc ; j++) \ |
|
5164
|
1272 { \ |
|
5275
|
1273 for (octave_idx_type i = 0 ; i < m1_nr ; i++) \ |
|
5164
|
1274 { \ |
|
|
1275 if ((m1.elem(i, j) OP m2.elem(i, j)) != ZERO) \ |
|
|
1276 { \ |
|
|
1277 r.data (ii) = m1.elem(i, j) OP m2.elem(i,j); \ |
|
|
1278 r.ridx (ii++) = i; \ |
|
|
1279 } \ |
|
|
1280 } \ |
|
|
1281 r.cidx(j+1) = ii; \ |
|
|
1282 } \ |
|
|
1283 } \ |
|
|
1284 \ |
|
|
1285 return r; \ |
|
|
1286 } |
|
|
1287 |
|
5775
|
1288 // FIXME Pass a specific ZERO value |
|
5164
|
1289 #define SPARSE_SMM_BIN_OPS(R1, R2, M1, M2) \ |
|
|
1290 SPARSE_SMM_BIN_OP_1 (R1, operator +, +, M1, M2) \ |
|
|
1291 SPARSE_SMM_BIN_OP_1 (R1, operator -, -, M1, M2) \ |
|
|
1292 SPARSE_SMM_BIN_OP_2 (R2, product, *, M1, M2, 0.0) \ |
|
|
1293 SPARSE_SMM_BIN_OP_2 (R2, quotient, /, M1, M2, 0.0) |
|
|
1294 |
|
6708
|
1295 #define SPARSE_SMM_CMP_OP_DECLS(M1, M2, API) \ |
|
|
1296 SPARSE_CMP_OP_DECL (mx_el_lt, M1, M2, API); \ |
|
|
1297 SPARSE_CMP_OP_DECL (mx_el_le, M1, M2, API); \ |
|
|
1298 SPARSE_CMP_OP_DECL (mx_el_ge, M1, M2, API); \ |
|
|
1299 SPARSE_CMP_OP_DECL (mx_el_gt, M1, M2, API); \ |
|
|
1300 SPARSE_CMP_OP_DECL (mx_el_eq, M1, M2, API); \ |
|
|
1301 SPARSE_CMP_OP_DECL (mx_el_ne, M1, M2, API); |
|
5164
|
1302 |
|
6708
|
1303 #define SPARSE_SMM_EQNE_OP_DECLS(M1, M2, API) \ |
|
|
1304 SPARSE_CMP_OP_DECL (mx_el_eq, M1, M2, API); \ |
|
|
1305 SPARSE_CMP_OP_DECL (mx_el_ne, M1, M2, API); |
|
5164
|
1306 |
|
|
1307 #define SPARSE_SMM_CMP_OP(F, OP, M1, C1, M2, C2) \ |
|
|
1308 SparseBoolMatrix \ |
|
|
1309 F (const M1& m1, const M2& m2) \ |
|
|
1310 { \ |
|
|
1311 SparseBoolMatrix r; \ |
|
|
1312 \ |
|
5275
|
1313 octave_idx_type m1_nr = m1.rows (); \ |
|
|
1314 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
1315 \ |
|
5275
|
1316 octave_idx_type m2_nr = m2.rows (); \ |
|
|
1317 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
1318 \ |
|
6221
|
1319 if (m1_nr == 1 && m1_nc == 1) \ |
|
|
1320 r = SparseBoolMatrix (F (m1.elem(0,0), m2)); \ |
|
|
1321 else if (m1_nr == m2_nr && m1_nc == m2_nc) \ |
|
5164
|
1322 { \ |
|
|
1323 if (m1_nr != 0 || m1_nc != 0) \ |
|
|
1324 { \ |
|
|
1325 /* Count num of non-zero elements */ \ |
|
5275
|
1326 octave_idx_type nel = 0; \ |
|
|
1327 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
|
1328 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1329 if (C1 (m1.elem(i, j)) OP C2 (m2.elem(i, j))) \ |
|
|
1330 nel++; \ |
|
|
1331 \ |
|
|
1332 r = SparseBoolMatrix (m1_nr, m1_nc, nel); \ |
|
|
1333 \ |
|
5275
|
1334 octave_idx_type ii = 0; \ |
|
5164
|
1335 r.cidx (0) = 0; \ |
|
5275
|
1336 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
5164
|
1337 { \ |
|
5275
|
1338 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1339 { \ |
|
|
1340 bool el = C1 (m1.elem(i, j)) OP C2 (m2.elem(i, j)); \ |
|
|
1341 if (el) \ |
|
|
1342 { \ |
|
|
1343 r.data(ii) = el; \ |
|
|
1344 r.ridx(ii++) = i; \ |
|
|
1345 } \ |
|
|
1346 } \ |
|
|
1347 r.cidx(j+1) = ii; \ |
|
|
1348 } \ |
|
|
1349 } \ |
|
|
1350 } \ |
|
|
1351 else \ |
|
|
1352 { \ |
|
|
1353 if ((m1_nr != 0 || m1_nc != 0) && (m2_nr != 0 || m2_nc != 0)) \ |
|
|
1354 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
1355 } \ |
|
|
1356 return r; \ |
|
|
1357 } |
|
|
1358 |
|
|
1359 #define SPARSE_SMM_CMP_OPS(M1, Z1, C1, M2, Z2, C2) \ |
|
|
1360 SPARSE_SMM_CMP_OP (mx_el_lt, <, M1, C1, M2, C2) \ |
|
|
1361 SPARSE_SMM_CMP_OP (mx_el_le, <=, M1, C1, M2, C2) \ |
|
|
1362 SPARSE_SMM_CMP_OP (mx_el_ge, >=, M1, C1, M2, C2) \ |
|
|
1363 SPARSE_SMM_CMP_OP (mx_el_gt, >, M1, C1, M2, C2) \ |
|
|
1364 SPARSE_SMM_CMP_OP (mx_el_eq, ==, M1, , M2, ) \ |
|
|
1365 SPARSE_SMM_CMP_OP (mx_el_ne, !=, M1, , M2, ) |
|
|
1366 |
|
|
1367 #define SPARSE_SMM_EQNE_OPS(M1, Z1, C1, M2, Z2, C2) \ |
|
|
1368 SPARSE_SMM_CMP_OP (mx_el_eq, ==, M1, , M2, ) \ |
|
|
1369 SPARSE_SMM_CMP_OP (mx_el_ne, !=, M1, , M2, ) |
|
|
1370 |
|
6708
|
1371 #define SPARSE_SMM_BOOL_OP_DECLS(M1, M2, API) \ |
|
|
1372 SPARSE_BOOL_OP_DECL (mx_el_and, M1, M2, API); \ |
|
|
1373 SPARSE_BOOL_OP_DECL (mx_el_or, M1, M2, API); |
|
5164
|
1374 |
|
|
1375 #define SPARSE_SMM_BOOL_OP(F, OP, M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
1376 SparseBoolMatrix \ |
|
|
1377 F (const M1& m1, const M2& m2) \ |
|
|
1378 { \ |
|
|
1379 SparseBoolMatrix r; \ |
|
|
1380 \ |
|
5275
|
1381 octave_idx_type m1_nr = m1.rows (); \ |
|
|
1382 octave_idx_type m1_nc = m1.cols (); \ |
|
5164
|
1383 \ |
|
5275
|
1384 octave_idx_type m2_nr = m2.rows (); \ |
|
|
1385 octave_idx_type m2_nc = m2.cols (); \ |
|
5164
|
1386 \ |
|
6221
|
1387 if (m1_nr == 1 && m1_nc == 1) \ |
|
|
1388 r = SparseBoolMatrix (F (m1.elem(0,0), m2)); \ |
|
|
1389 else if (m1_nr == m2_nr && m1_nc == m2_nc) \ |
|
5164
|
1390 { \ |
|
|
1391 if (m1_nr != 0 || m1_nc != 0) \ |
|
|
1392 { \ |
|
|
1393 /* Count num of non-zero elements */ \ |
|
5275
|
1394 octave_idx_type nel = 0; \ |
|
|
1395 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
|
1396 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1397 if ((m1.elem(i, j) != LHS_ZERO) \ |
|
|
1398 OP (m2.elem(i, j) != RHS_ZERO)) \ |
|
|
1399 nel++; \ |
|
|
1400 \ |
|
|
1401 r = SparseBoolMatrix (m1_nr, m1_nc, nel); \ |
|
|
1402 \ |
|
5275
|
1403 octave_idx_type ii = 0; \ |
|
5164
|
1404 r.cidx (0) = 0; \ |
|
5275
|
1405 for (octave_idx_type j = 0; j < m1_nc; j++) \ |
|
5164
|
1406 { \ |
|
5275
|
1407 for (octave_idx_type i = 0; i < m1_nr; i++) \ |
|
5164
|
1408 { \ |
|
|
1409 bool el = (m1.elem(i, j) != LHS_ZERO) \ |
|
|
1410 OP (m2.elem(i, j) != RHS_ZERO); \ |
|
|
1411 if (el) \ |
|
|
1412 { \ |
|
|
1413 r.data(ii) = el; \ |
|
|
1414 r.ridx(ii++) = i; \ |
|
|
1415 } \ |
|
|
1416 } \ |
|
|
1417 r.cidx(j+1) = ii; \ |
|
|
1418 } \ |
|
|
1419 } \ |
|
|
1420 } \ |
|
|
1421 else \ |
|
|
1422 { \ |
|
|
1423 if ((m1_nr != 0 || m1_nc != 0) && (m2_nr != 0 || m2_nc != 0)) \ |
|
|
1424 gripe_nonconformant (#F, m1_nr, m1_nc, m2_nr, m2_nc); \ |
|
|
1425 } \ |
|
|
1426 return r; \ |
|
|
1427 } |
|
|
1428 |
|
|
1429 #define SPARSE_SMM_BOOL_OPS2(M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
1430 SPARSE_SMM_BOOL_OP (mx_el_and, &&, M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
1431 SPARSE_SMM_BOOL_OP (mx_el_or, ||, M1, M2, LHS_ZERO, RHS_ZERO) \ |
|
|
1432 |
|
|
1433 #define SPARSE_SMM_BOOL_OPS(M1, M2, ZERO) \ |
|
|
1434 SPARSE_SMM_BOOL_OPS2(M1, M2, ZERO, ZERO) |
|
|
1435 |
|
6708
|
1436 #define SPARSE_SMM_OP_DECLS(R1, R2, M1, M2, API) \ |
|
|
1437 SPARSE_SMM_BIN_OP_DECLS (R1, R2, M1, M2, API) \ |
|
|
1438 SPARSE_SMM_CMP_OP_DECLS (M1, M2, API) \ |
|
|
1439 SPARSE_SMM_BOOL_OP_DECLS (M1, M2, API) |
|
5164
|
1440 |
|
|
1441 // Avoid some code duplication. Maybe we should use templates. |
|
|
1442 |
|
|
1443 #define SPARSE_CUMSUM(RET_TYPE, ELT_TYPE, FCN) \ |
|
|
1444 \ |
|
5275
|
1445 octave_idx_type nr = rows (); \ |
|
|
1446 octave_idx_type nc = cols (); \ |
|
5164
|
1447 \ |
|
|
1448 RET_TYPE retval; \ |
|
|
1449 \ |
|
|
1450 if (nr > 0 && nc > 0) \ |
|
|
1451 { \ |
|
|
1452 if ((nr == 1 && dim == -1) || dim == 1) \ |
|
|
1453 /* Ugly!! Is there a better way? */ \ |
|
|
1454 retval = transpose (). FCN (0) .transpose (); \ |
|
|
1455 else \ |
|
|
1456 { \ |
|
5275
|
1457 octave_idx_type nel = 0; \ |
|
|
1458 for (octave_idx_type i = 0; i < nc; i++) \ |
|
5164
|
1459 { \ |
|
|
1460 ELT_TYPE t = ELT_TYPE (); \ |
|
5275
|
1461 for (octave_idx_type j = cidx (i); j < cidx (i+1); j++) \ |
|
5164
|
1462 { \ |
|
|
1463 t += data(j); \ |
|
|
1464 if (t != ELT_TYPE ()) \ |
|
6482
|
1465 { \ |
|
|
1466 if (j == cidx(i+1) - 1) \ |
|
|
1467 nel += nr - ridx(j); \ |
|
|
1468 else \ |
|
|
1469 nel += ridx(j+1) - ridx(j); \ |
|
|
1470 } \ |
|
5164
|
1471 } \ |
|
|
1472 } \ |
|
|
1473 retval = RET_TYPE (nr, nc, nel); \ |
|
|
1474 retval.cidx(0) = 0; \ |
|
5275
|
1475 octave_idx_type ii = 0; \ |
|
|
1476 for (octave_idx_type i = 0; i < nc; i++) \ |
|
5164
|
1477 { \ |
|
|
1478 ELT_TYPE t = ELT_TYPE (); \ |
|
5275
|
1479 for (octave_idx_type j = cidx (i); j < cidx (i+1); j++) \ |
|
5164
|
1480 { \ |
|
|
1481 t += data(j); \ |
|
|
1482 if (t != ELT_TYPE ()) \ |
|
|
1483 { \ |
|
|
1484 if (j == cidx(i+1) - 1) \ |
|
|
1485 { \ |
|
5275
|
1486 for (octave_idx_type k = ridx(j); k < nr; k++) \ |
|
5164
|
1487 { \ |
|
|
1488 retval.data (ii) = t; \ |
|
|
1489 retval.ridx (ii++) = k; \ |
|
|
1490 } \ |
|
|
1491 } \ |
|
|
1492 else \ |
|
|
1493 { \ |
|
5275
|
1494 for (octave_idx_type k = ridx(j); k < ridx(j+1); k++) \ |
|
5164
|
1495 { \ |
|
|
1496 retval.data (ii) = t; \ |
|
|
1497 retval.ridx (ii++) = k; \ |
|
|
1498 } \ |
|
|
1499 } \ |
|
|
1500 } \ |
|
|
1501 } \ |
|
|
1502 retval.cidx(i+1) = ii; \ |
|
|
1503 } \ |
|
|
1504 } \ |
|
|
1505 } \ |
|
|
1506 else \ |
|
|
1507 retval = RET_TYPE (nr,nc); \ |
|
|
1508 \ |
|
|
1509 return retval |
|
|
1510 |
|
|
1511 |
|
|
1512 #define SPARSE_CUMPROD(RET_TYPE, ELT_TYPE, FCN) \ |
|
|
1513 \ |
|
5275
|
1514 octave_idx_type nr = rows (); \ |
|
|
1515 octave_idx_type nc = cols (); \ |
|
5164
|
1516 \ |
|
|
1517 RET_TYPE retval; \ |
|
|
1518 \ |
|
|
1519 if (nr > 0 && nc > 0) \ |
|
|
1520 { \ |
|
|
1521 if ((nr == 1 && dim == -1) || dim == 1) \ |
|
|
1522 /* Ugly!! Is there a better way? */ \ |
|
|
1523 retval = transpose (). FCN (0) .transpose (); \ |
|
|
1524 else \ |
|
|
1525 { \ |
|
5275
|
1526 octave_idx_type nel = 0; \ |
|
|
1527 for (octave_idx_type i = 0; i < nc; i++) \ |
|
5164
|
1528 { \ |
|
5275
|
1529 octave_idx_type jj = 0; \ |
|
|
1530 for (octave_idx_type j = cidx (i); j < cidx (i+1); j++) \ |
|
5164
|
1531 { \ |
|
|
1532 if (jj == ridx(j)) \ |
|
|
1533 { \ |
|
|
1534 nel++; \ |
|
|
1535 jj++; \ |
|
|
1536 } \ |
|
|
1537 else \ |
|
|
1538 break; \ |
|
|
1539 } \ |
|
|
1540 } \ |
|
|
1541 retval = RET_TYPE (nr, nc, nel); \ |
|
|
1542 retval.cidx(0) = 0; \ |
|
5275
|
1543 octave_idx_type ii = 0; \ |
|
|
1544 for (octave_idx_type i = 0; i < nc; i++) \ |
|
5164
|
1545 { \ |
|
|
1546 ELT_TYPE t = ELT_TYPE (1.); \ |
|
5275
|
1547 octave_idx_type jj = 0; \ |
|
|
1548 for (octave_idx_type j = cidx (i); j < cidx (i+1); j++) \ |
|
5164
|
1549 { \ |
|
|
1550 if (jj == ridx(j)) \ |
|
|
1551 { \ |
|
|
1552 t *= data(j); \ |
|
|
1553 retval.data(ii) = t; \ |
|
|
1554 retval.ridx(ii++) = jj++; \ |
|
|
1555 } \ |
|
|
1556 else \ |
|
|
1557 break; \ |
|
|
1558 } \ |
|
|
1559 retval.cidx(i+1) = ii; \ |
|
|
1560 } \ |
|
|
1561 } \ |
|
|
1562 } \ |
|
|
1563 else \ |
|
|
1564 retval = RET_TYPE (nr,nc); \ |
|
|
1565 \ |
|
|
1566 return retval |
|
|
1567 |
|
|
1568 #define SPARSE_BASE_REDUCTION_OP(RET_TYPE, EL_TYPE, ROW_EXPR, COL_EXPR, \ |
|
|
1569 INIT_VAL, MT_RESULT) \ |
|
|
1570 \ |
|
5275
|
1571 octave_idx_type nr = rows (); \ |
|
|
1572 octave_idx_type nc = cols (); \ |
|
5164
|
1573 \ |
|
|
1574 RET_TYPE retval; \ |
|
|
1575 \ |
|
|
1576 if (nr > 0 && nc > 0) \ |
|
|
1577 { \ |
|
|
1578 if ((nr == 1 && dim == -1) || dim == 1) \ |
|
|
1579 { \ |
|
|
1580 OCTAVE_LOCAL_BUFFER (EL_TYPE, tmp, nr); \ |
|
|
1581 \ |
|
5275
|
1582 for (octave_idx_type i = 0; i < nr; i++) \ |
|
5164
|
1583 { \ |
|
|
1584 tmp[i] = INIT_VAL; \ |
|
5275
|
1585 for (octave_idx_type j = 0; j < nc; j++) \ |
|
5164
|
1586 { \ |
|
|
1587 ROW_EXPR; \ |
|
|
1588 } \ |
|
|
1589 } \ |
|
5275
|
1590 octave_idx_type nel = 0; \ |
|
|
1591 for (octave_idx_type i = 0; i < nr; i++) \ |
|
5164
|
1592 if (tmp[i] != EL_TYPE ()) \ |
|
|
1593 nel++ ; \ |
|
5275
|
1594 retval = RET_TYPE (nr, static_cast<octave_idx_type> (1), nel); \ |
|
5164
|
1595 retval.cidx(0) = 0; \ |
|
|
1596 retval.cidx(1) = nel; \ |
|
|
1597 nel = 0; \ |
|
5275
|
1598 for (octave_idx_type i = 0; i < nr; i++) \ |
|
5164
|
1599 if (tmp[i] != EL_TYPE ()) \ |
|
|
1600 { \ |
|
|
1601 retval.data(nel) = tmp[i]; \ |
|
|
1602 retval.ridx(nel++) = i; \ |
|
|
1603 } \ |
|
|
1604 } \ |
|
|
1605 else \ |
|
|
1606 { \ |
|
|
1607 OCTAVE_LOCAL_BUFFER (EL_TYPE, tmp, nc); \ |
|
|
1608 \ |
|
5275
|
1609 for (octave_idx_type j = 0; j < nc; j++) \ |
|
5164
|
1610 { \ |
|
|
1611 tmp[j] = INIT_VAL; \ |
|
5275
|
1612 for (octave_idx_type i = 0; i < nr; i++) \ |
|
5164
|
1613 { \ |
|
|
1614 COL_EXPR; \ |
|
|
1615 } \ |
|
|
1616 } \ |
|
5275
|
1617 octave_idx_type nel = 0; \ |
|
|
1618 for (octave_idx_type i = 0; i < nc; i++) \ |
|
5164
|
1619 if (tmp[i] != EL_TYPE ()) \ |
|
|
1620 nel++ ; \ |
|
5275
|
1621 retval = RET_TYPE (static_cast<octave_idx_type> (1), nc, nel); \ |
|
5164
|
1622 retval.cidx(0) = 0; \ |
|
|
1623 nel = 0; \ |
|
5275
|
1624 for (octave_idx_type i = 0; i < nc; i++) \ |
|
5164
|
1625 if (tmp[i] != EL_TYPE ()) \ |
|
|
1626 { \ |
|
|
1627 retval.data(nel) = tmp[i]; \ |
|
|
1628 retval.ridx(nel++) = 0; \ |
|
|
1629 retval.cidx(i+1) = retval.cidx(i) + 1; \ |
|
|
1630 } \ |
|
|
1631 else \ |
|
|
1632 retval.cidx(i+1) = retval.cidx(i); \ |
|
|
1633 } \ |
|
|
1634 } \ |
|
|
1635 else if (nc == 0 && (nr == 0 || (nr == 1 && dim == -1))) \ |
|
|
1636 { \ |
|
5275
|
1637 retval = RET_TYPE (static_cast<octave_idx_type> (1), \ |
|
|
1638 static_cast<octave_idx_type> (1), \ |
|
|
1639 static_cast<octave_idx_type> (1)); \ |
|
5164
|
1640 retval.cidx(0) = 0; \ |
|
|
1641 retval.cidx(1) = 1; \ |
|
|
1642 retval.ridx(0) = 0; \ |
|
|
1643 retval.data(0) = MT_RESULT; \ |
|
|
1644 } \ |
|
|
1645 else if (nr == 0 && (dim == 0 || dim == -1)) \ |
|
|
1646 { \ |
|
5275
|
1647 retval = RET_TYPE (static_cast<octave_idx_type> (1), nc, nc); \ |
|
5164
|
1648 retval.cidx (0) = 0; \ |
|
5275
|
1649 for (octave_idx_type i = 0; i < nc ; i++) \ |
|
5164
|
1650 { \ |
|
|
1651 retval.ridx (i) = 0; \ |
|
|
1652 retval.cidx (i+1) = i; \ |
|
|
1653 retval.data (i) = MT_RESULT; \ |
|
|
1654 } \ |
|
|
1655 } \ |
|
|
1656 else if (nc == 0 && dim == 1) \ |
|
|
1657 { \ |
|
5275
|
1658 retval = RET_TYPE (nr, static_cast<octave_idx_type> (1), nr); \ |
|
5164
|
1659 retval.cidx(0) = 0; \ |
|
|
1660 retval.cidx(1) = nr; \ |
|
5275
|
1661 for (octave_idx_type i = 0; i < nr; i++) \ |
|
5164
|
1662 { \ |
|
|
1663 retval.ridx(i) = i; \ |
|
|
1664 retval.data(i) = MT_RESULT; \ |
|
|
1665 } \ |
|
|
1666 } \ |
|
|
1667 else \ |
|
|
1668 retval.resize (nr > 0, nc > 0); \ |
|
|
1669 \ |
|
|
1670 return retval |
|
|
1671 |
|
|
1672 #define SPARSE_REDUCTION_OP_ROW_EXPR(OP) \ |
|
|
1673 tmp[i] OP elem (i, j) |
|
|
1674 |
|
|
1675 #define SPARSE_REDUCTION_OP_COL_EXPR(OP) \ |
|
|
1676 tmp[j] OP elem (i, j) |
|
|
1677 |
|
|
1678 #define SPARSE_REDUCTION_OP(RET_TYPE, EL_TYPE, OP, INIT_VAL, MT_RESULT) \ |
|
|
1679 SPARSE_BASE_REDUCTION_OP (RET_TYPE, EL_TYPE, \ |
|
|
1680 SPARSE_REDUCTION_OP_ROW_EXPR (OP), \ |
|
|
1681 SPARSE_REDUCTION_OP_COL_EXPR (OP), \ |
|
|
1682 INIT_VAL, MT_RESULT) |
|
|
1683 |
|
|
1684 #define SPARSE_ANY_ALL_OP_ROW_CODE(TEST_OP, TEST_TRUE_VAL) \ |
|
|
1685 if (elem (i, j) TEST_OP 0.0) \ |
|
|
1686 { \ |
|
|
1687 tmp[i] = TEST_TRUE_VAL; \ |
|
|
1688 break; \ |
|
|
1689 } |
|
|
1690 |
|
|
1691 #define SPARSE_ANY_ALL_OP_COL_CODE(TEST_OP, TEST_TRUE_VAL) \ |
|
|
1692 if (elem (i, j) TEST_OP 0.0) \ |
|
|
1693 { \ |
|
|
1694 tmp[j] = TEST_TRUE_VAL; \ |
|
|
1695 break; \ |
|
|
1696 } |
|
|
1697 |
|
|
1698 #define SPARSE_ANY_ALL_OP(DIM, INIT_VAL, TEST_OP, TEST_TRUE_VAL) \ |
|
|
1699 SPARSE_BASE_REDUCTION_OP (SparseBoolMatrix, char, \ |
|
|
1700 SPARSE_ANY_ALL_OP_ROW_CODE (TEST_OP, TEST_TRUE_VAL), \ |
|
|
1701 SPARSE_ANY_ALL_OP_COL_CODE (TEST_OP, TEST_TRUE_VAL), \ |
|
|
1702 INIT_VAL, INIT_VAL) |
|
|
1703 |
|
|
1704 #define SPARSE_ALL_OP(DIM) SPARSE_ANY_ALL_OP (DIM, true, ==, false) |
|
|
1705 |
|
|
1706 #define SPARSE_ANY_OP(DIM) SPARSE_ANY_ALL_OP (DIM, false, !=, true) |
|
|
1707 |
|
5681
|
1708 #define SPARSE_SPARSE_MUL( RET_TYPE, RET_EL_TYPE, EL_TYPE ) \ |
|
5275
|
1709 octave_idx_type nr = m.rows (); \ |
|
|
1710 octave_idx_type nc = m.cols (); \ |
|
5164
|
1711 \ |
|
5275
|
1712 octave_idx_type a_nr = a.rows (); \ |
|
|
1713 octave_idx_type a_nc = a.cols (); \ |
|
5164
|
1714 \ |
|
6221
|
1715 if (nr == 1 && nc == 1) \ |
|
|
1716 { \ |
|
|
1717 RET_EL_TYPE s = m.elem(0,0); \ |
|
|
1718 octave_idx_type nz = a.nnz(); \ |
|
|
1719 RET_TYPE r (a_nr, a_nc, nz); \ |
|
|
1720 \ |
|
|
1721 for (octave_idx_type i = 0; i < nz; i++) \ |
|
|
1722 { \ |
|
|
1723 OCTAVE_QUIT; \ |
|
|
1724 r.data(i) = s * a.data(i); \ |
|
|
1725 r.ridx(i) = a.ridx(i); \ |
|
|
1726 } \ |
|
|
1727 for (octave_idx_type i = 0; i < a_nc + 1; i++) \ |
|
|
1728 { \ |
|
|
1729 OCTAVE_QUIT; \ |
|
|
1730 r.cidx(i) = a.cidx(i); \ |
|
|
1731 } \ |
|
|
1732 \ |
|
|
1733 r.maybe_compress (true); \ |
|
|
1734 return r; \ |
|
|
1735 } \ |
|
|
1736 else if (a_nr == 1 && a_nc == 1) \ |
|
|
1737 { \ |
|
|
1738 RET_EL_TYPE s = a.elem(0,0); \ |
|
|
1739 octave_idx_type nz = m.nnz(); \ |
|
|
1740 RET_TYPE r (nr, nc, nz); \ |
|
|
1741 \ |
|
|
1742 for (octave_idx_type i = 0; i < nz; i++) \ |
|
|
1743 { \ |
|
|
1744 OCTAVE_QUIT; \ |
|
|
1745 r.data(i) = m.data(i) * s; \ |
|
|
1746 r.ridx(i) = m.ridx(i); \ |
|
|
1747 } \ |
|
|
1748 for (octave_idx_type i = 0; i < nc + 1; i++) \ |
|
|
1749 { \ |
|
|
1750 OCTAVE_QUIT; \ |
|
|
1751 r.cidx(i) = m.cidx(i); \ |
|
|
1752 } \ |
|
|
1753 \ |
|
|
1754 r.maybe_compress (true); \ |
|
|
1755 return r; \ |
|
|
1756 } \ |
|
|
1757 else if (nc != a_nr) \ |
|
5164
|
1758 { \ |
|
|
1759 gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); \ |
|
|
1760 return RET_TYPE (); \ |
|
|
1761 } \ |
|
|
1762 else \ |
|
|
1763 { \ |
|
5586
|
1764 OCTAVE_LOCAL_BUFFER (octave_idx_type, w, nr); \ |
|
5876
|
1765 RET_TYPE retval (nr, a_nc, static_cast<octave_idx_type> (0)); \ |
|
5586
|
1766 for (octave_idx_type i = 0; i < nr; i++) \ |
|
|
1767 w[i] = 0; \ |
|
5795
|
1768 retval.xcidx(0) = 0; \ |
|
5164
|
1769 \ |
|
5275
|
1770 octave_idx_type nel = 0; \ |
|
5164
|
1771 \ |
|
5275
|
1772 for (octave_idx_type i = 0; i < a_nc; i++) \ |
|
5164
|
1773 { \ |
|
5275
|
1774 for (octave_idx_type j = a.cidx(i); j < a.cidx(i+1); j++) \ |
|
5164
|
1775 { \ |
|
5275
|
1776 octave_idx_type col = a.ridx(j); \ |
|
|
1777 for (octave_idx_type k = m.cidx(col) ; k < m.cidx(col+1); k++) \ |
|
5586
|
1778 { \ |
|
|
1779 if (w[m.ridx(k)] < i + 1) \ |
|
|
1780 { \ |
|
|
1781 w[m.ridx(k)] = i + 1; \ |
|
|
1782 nel++; \ |
|
|
1783 } \ |
|
5587
|
1784 OCTAVE_QUIT; \ |
|
5586
|
1785 } \ |
|
5164
|
1786 } \ |
|
5795
|
1787 retval.xcidx(i+1) = nel; \ |
|
5164
|
1788 } \ |
|
|
1789 \ |
|
|
1790 if (nel == 0) \ |
|
|
1791 return RET_TYPE (nr, a_nc); \ |
|
|
1792 else \ |
|
|
1793 { \ |
|
5586
|
1794 for (octave_idx_type i = 0; i < nr; i++) \ |
|
|
1795 w[i] = 0; \ |
|
|
1796 \ |
|
5681
|
1797 OCTAVE_LOCAL_BUFFER (RET_EL_TYPE, Xcol, nr); \ |
|
5586
|
1798 \ |
|
5795
|
1799 retval.change_capacity (nel); \ |
|
5587
|
1800 /* The optimal break-point as estimated from simulations */ \ |
|
|
1801 /* Note that Mergesort is O(nz log(nz)) while searching all */ \ |
|
|
1802 /* values is O(nr), where nz here is non-zero per row of */ \ |
|
|
1803 /* length nr. The test itself was then derived from the */ \ |
|
|
1804 /* simulation with random square matrices and the observation */ \ |
|
|
1805 /* of the number of non-zero elements in the output matrix */ \ |
|
|
1806 /* it was found that the breakpoints were */ \ |
|
|
1807 /* nr: 500 1000 2000 5000 10000 */ \ |
|
|
1808 /* nz: 6 25 97 585 2202 */ \ |
|
|
1809 /* The below is a simplication of the 'polyfit'-ed parameters */ \ |
|
|
1810 /* to these breakpoints */ \ |
|
5795
|
1811 octave_idx_type n_per_col = (a_nc > 43000 ? 43000 : \ |
|
|
1812 (a_nc * a_nc) / 43000); \ |
|
|
1813 octave_idx_type ii = 0; \ |
|
|
1814 octave_idx_type *ri = retval.xridx(); \ |
|
|
1815 octave_sort<octave_idx_type> sort; \ |
|
|
1816 \ |
|
|
1817 for (octave_idx_type i = 0; i < a_nc ; i++) \ |
|
5164
|
1818 { \ |
|
5795
|
1819 if (retval.xcidx(i+1) - retval.xcidx(i) > n_per_col) \ |
|
5587
|
1820 { \ |
|
|
1821 for (octave_idx_type j = a.cidx(i); j < a.cidx(i+1); j++) \ |
|
|
1822 { \ |
|
|
1823 octave_idx_type col = a.ridx(j); \ |
|
|
1824 EL_TYPE tmpval = a.data(j); \ |
|
|
1825 for (octave_idx_type k = m.cidx(col) ; \ |
|
|
1826 k < m.cidx(col+1); k++) \ |
|
|
1827 { \ |
|
|
1828 OCTAVE_QUIT; \ |
|
|
1829 octave_idx_type row = m.ridx(k); \ |
|
|
1830 if (w[row] < i + 1) \ |
|
|
1831 { \ |
|
|
1832 w[row] = i + 1; \ |
|
|
1833 Xcol[row] = tmpval * m.data(k); \ |
|
|
1834 } \ |
|
|
1835 else \ |
|
|
1836 Xcol[row] += tmpval * m.data(k); \ |
|
|
1837 } \ |
|
|
1838 } \ |
|
|
1839 for (octave_idx_type k = 0; k < nr; k++) \ |
|
5813
|
1840 if (w[k] == i + 1) \ |
|
5587
|
1841 { \ |
|
|
1842 retval.xdata(ii) = Xcol[k]; \ |
|
|
1843 retval.xridx(ii++) = k; \ |
|
|
1844 } \ |
|
5795
|
1845 } \ |
|
|
1846 else \ |
|
|
1847 { \ |
|
|
1848 for (octave_idx_type j = a.cidx(i); j < a.cidx(i+1); j++) \ |
|
|
1849 { \ |
|
|
1850 octave_idx_type col = a.ridx(j); \ |
|
|
1851 EL_TYPE tmpval = a.data(j); \ |
|
|
1852 for (octave_idx_type k = m.cidx(col) ; \ |
|
|
1853 k < m.cidx(col+1); k++) \ |
|
|
1854 { \ |
|
|
1855 OCTAVE_QUIT; \ |
|
|
1856 octave_idx_type row = m.ridx(k); \ |
|
|
1857 if (w[row] < i + 1) \ |
|
|
1858 { \ |
|
|
1859 w[row] = i + 1; \ |
|
|
1860 retval.xridx(ii++) = row;\ |
|
|
1861 Xcol[row] = tmpval * m.data(k); \ |
|
|
1862 } \ |
|
|
1863 else \ |
|
|
1864 Xcol[row] += tmpval * m.data(k); \ |
|
|
1865 } \ |
|
|
1866 } \ |
|
|
1867 sort.sort (ri + retval.xcidx(i), ii - retval.xcidx(i)); \ |
|
|
1868 for (octave_idx_type k = retval.xcidx(i); k < ii; k++) \ |
|
|
1869 retval.xdata(k) = Xcol[retval.xridx(k)]; \ |
|
5587
|
1870 } \ |
|
5164
|
1871 } \ |
|
5813
|
1872 retval.maybe_compress (true);\ |
|
5164
|
1873 return retval; \ |
|
|
1874 } \ |
|
|
1875 } |
|
|
1876 |
|
5681
|
1877 #define SPARSE_FULL_MUL( RET_TYPE, EL_TYPE, ZERO ) \ |
|
5429
|
1878 octave_idx_type nr = m.rows (); \ |
|
|
1879 octave_idx_type nc = m.cols (); \ |
|
|
1880 \ |
|
|
1881 octave_idx_type a_nr = a.rows (); \ |
|
|
1882 octave_idx_type a_nc = a.cols (); \ |
|
|
1883 \ |
|
6221
|
1884 if (nr == 1 && nc == 1) \ |
|
|
1885 { \ |
|
|
1886 RET_TYPE retval (a_nr, a_nc, ZERO); \ |
|
|
1887 for (octave_idx_type i = 0; i < a_nc ; i++) \ |
|
|
1888 { \ |
|
|
1889 for (octave_idx_type j = 0; j < a_nr; j++) \ |
|
|
1890 { \ |
|
|
1891 OCTAVE_QUIT; \ |
|
|
1892 retval.elem (j,i) += a.elem(j,i) * m.elem(0,0); \ |
|
|
1893 } \ |
|
|
1894 } \ |
|
|
1895 return retval; \ |
|
|
1896 } \ |
|
|
1897 else if (nc != a_nr) \ |
|
5429
|
1898 { \ |
|
|
1899 gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); \ |
|
|
1900 return RET_TYPE (); \ |
|
|
1901 } \ |
|
|
1902 else \ |
|
|
1903 { \ |
|
5681
|
1904 RET_TYPE retval (nr, a_nc, ZERO); \ |
|
5429
|
1905 \ |
|
|
1906 for (octave_idx_type i = 0; i < a_nc ; i++) \ |
|
|
1907 { \ |
|
|
1908 for (octave_idx_type j = 0; j < a_nr; j++) \ |
|
|
1909 { \ |
|
|
1910 OCTAVE_QUIT; \ |
|
|
1911 \ |
|
|
1912 EL_TYPE tmpval = a.elem(j,i); \ |
|
|
1913 for (octave_idx_type k = m.cidx(j) ; k < m.cidx(j+1); k++) \ |
|
|
1914 retval.elem (m.ridx(k),i) += tmpval * m.data(k); \ |
|
|
1915 } \ |
|
|
1916 } \ |
|
|
1917 return retval; \ |
|
|
1918 } |
|
|
1919 |
|
5681
|
1920 #define FULL_SPARSE_MUL( RET_TYPE, EL_TYPE, ZERO ) \ |
|
5429
|
1921 octave_idx_type nr = m.rows (); \ |
|
|
1922 octave_idx_type nc = m.cols (); \ |
|
|
1923 \ |
|
|
1924 octave_idx_type a_nr = a.rows (); \ |
|
|
1925 octave_idx_type a_nc = a.cols (); \ |
|
|
1926 \ |
|
6221
|
1927 if (a_nr == 1 && a_nc == 1) \ |
|
|
1928 { \ |
|
|
1929 RET_TYPE retval (nr, nc, ZERO); \ |
|
|
1930 for (octave_idx_type i = 0; i < nc ; i++) \ |
|
|
1931 { \ |
|
|
1932 for (octave_idx_type j = 0; j < nr; j++) \ |
|
|
1933 { \ |
|
|
1934 OCTAVE_QUIT; \ |
|
|
1935 retval.elem (j,i) += a.elem(0,0) * m.elem(j,i); \ |
|
|
1936 } \ |
|
|
1937 } \ |
|
|
1938 return retval; \ |
|
|
1939 } \ |
|
|
1940 else if (nc != a_nr) \ |
|
5429
|
1941 { \ |
|
|
1942 gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); \ |
|
|
1943 return RET_TYPE (); \ |
|
|
1944 } \ |
|
|
1945 else \ |
|
|
1946 { \ |
|
5681
|
1947 RET_TYPE retval (nr, a_nc, ZERO); \ |
|
5429
|
1948 \ |
|
|
1949 for (octave_idx_type i = 0; i < a_nc ; i++) \ |
|
|
1950 { \ |
|
|
1951 for (octave_idx_type j = a.cidx(i); j < a.cidx(i+1); j++) \ |
|
|
1952 { \ |
|
|
1953 octave_idx_type col = a.ridx(j); \ |
|
|
1954 EL_TYPE tmpval = a.data(j); \ |
|
|
1955 OCTAVE_QUIT; \ |
|
|
1956 \ |
|
|
1957 for (octave_idx_type k = 0 ; k < nr; k++) \ |
|
|
1958 retval.elem (k,i) += tmpval * m.elem(k,col); \ |
|
|
1959 } \ |
|
|
1960 } \ |
|
|
1961 return retval; \ |
|
|
1962 } |
|
|
1963 |
|
5164
|
1964 #endif |
|
|
1965 |
|
|
1966 /* |
|
|
1967 ;;; Local Variables: *** |
|
|
1968 ;;; mode: C++ *** |
|
|
1969 ;;; End: *** |
|
|
1970 */ |
