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