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