Mercurial > hg > octave-lyh

comparison liboctave/MSparse.cc @ 5275:23b37da9fd5b

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 2005-04-08 16:07:35 by jwe]
author jwe
date Fri, 08 Apr 2005 16:07:37 +0000
parents 57077d0ddc8e
children 4c8a2e4e0717
comparison
equal deleted inserted replaced
5274:eae7b40388e9 5275:23b37da9fd5b
39 MSparse<T>& 39 MSparse<T>&
40 operator += (MSparse<T>& a, const MSparse<T>& b) 40 operator += (MSparse<T>& a, const MSparse<T>& b)
41 { 41 {
42 MSparse<T> r; 42 MSparse<T> r;
43 43
44 int a_nr = a.rows (); 44 octave_idx_type a_nr = a.rows ();
45 int a_nc = a.cols (); 45 octave_idx_type a_nc = a.cols ();
46 46
47 int b_nr = b.rows (); 47 octave_idx_type b_nr = b.rows ();
48 int b_nc = b.cols (); 48 octave_idx_type b_nc = b.cols ();
49 49
50 if (a_nr != b_nr || a_nc != b_nc) 50 if (a_nr != b_nr || a_nc != b_nc)
51 gripe_nonconformant ("operator +=" , a_nr, a_nc, b_nr, b_nc); 51 gripe_nonconformant ("operator +=" , a_nr, a_nc, b_nr, b_nc);
52 else 52 else
53 { 53 {
54 r = MSparse<T> (a_nr, a_nc, (a.nnz () + b.nnz ())); 54 r = MSparse<T> (a_nr, a_nc, (a.nnz () + b.nnz ()));
55 55
56 int jx = 0; 56 octave_idx_type jx = 0;
57 for (int i = 0 ; i < a_nc ; i++) 57 for (octave_idx_type i = 0 ; i < a_nc ; i++)
58 { 58 {
59 int ja = a.cidx(i); 59 octave_idx_type ja = a.cidx(i);
60 int ja_max = a.cidx(i+1); 60 octave_idx_type ja_max = a.cidx(i+1);
61 bool ja_lt_max= ja < ja_max; 61 bool ja_lt_max= ja < ja_max;
62 62
63 int jb = b.cidx(i); 63 octave_idx_type jb = b.cidx(i);
64 int jb_max = b.cidx(i+1); 64 octave_idx_type jb_max = b.cidx(i+1);
65 bool jb_lt_max = jb < jb_max; 65 bool jb_lt_max = jb < jb_max;
66 66
67 while (ja_lt_max || jb_lt_max ) 67 while (ja_lt_max || jb_lt_max )
68 { 68 {
69 OCTAVE_QUIT; 69 OCTAVE_QUIT;
112 MSparse<T>& 112 MSparse<T>&
113 operator -= (MSparse<T>& a, const MSparse<T>& b) 113 operator -= (MSparse<T>& a, const MSparse<T>& b)
114 { 114 {
115 MSparse<T> r; 115 MSparse<T> r;
116 116
117 int a_nr = a.rows (); 117 octave_idx_type a_nr = a.rows ();
118 int a_nc = a.cols (); 118 octave_idx_type a_nc = a.cols ();
119 119
120 int b_nr = b.rows (); 120 octave_idx_type b_nr = b.rows ();
121 int b_nc = b.cols (); 121 octave_idx_type b_nc = b.cols ();
122 122
123 if (a_nr != b_nr || a_nc != b_nc) 123 if (a_nr != b_nr || a_nc != b_nc)
124 gripe_nonconformant ("operator -=" , a_nr, a_nc, b_nr, b_nc); 124 gripe_nonconformant ("operator -=" , a_nr, a_nc, b_nr, b_nc);
125 else 125 else
126 { 126 {
127 r = MSparse<T> (a_nr, a_nc, (a.nnz () + b.nnz ())); 127 r = MSparse<T> (a_nr, a_nc, (a.nnz () + b.nnz ()));
128 128
129 int jx = 0; 129 octave_idx_type jx = 0;
130 for (int i = 0 ; i < a_nc ; i++) 130 for (octave_idx_type i = 0 ; i < a_nc ; i++)
131 { 131 {
132 int ja = a.cidx(i); 132 octave_idx_type ja = a.cidx(i);
133 int ja_max = a.cidx(i+1); 133 octave_idx_type ja_max = a.cidx(i+1);
134 bool ja_lt_max= ja < ja_max; 134 bool ja_lt_max= ja < ja_max;
135 135
136 int jb = b.cidx(i); 136 octave_idx_type jb = b.cidx(i);
137 int jb_max = b.cidx(i+1); 137 octave_idx_type jb_max = b.cidx(i+1);
138 bool jb_lt_max = jb < jb_max; 138 bool jb_lt_max = jb < jb_max;
139 139
140 while (ja_lt_max || jb_lt_max ) 140 while (ja_lt_max || jb_lt_max )
141 { 141 {
142 OCTAVE_QUIT; 142 OCTAVE_QUIT;
186 #define SPARSE_A2S_OP_1(OP) \ 186 #define SPARSE_A2S_OP_1(OP) \
187 template <class T> \ 187 template <class T> \
188 MArray2<T> \ 188 MArray2<T> \
189 operator OP (const MSparse<T>& a, const T& s) \ 189 operator OP (const MSparse<T>& a, const T& s) \
190 { \ 190 { \
191 int nr = a.rows (); \ 191 octave_idx_type nr = a.rows (); \
192 int nc = a.cols (); \ 192 octave_idx_type nc = a.cols (); \
193 \ 193 \
194 MArray2<T> r (nr, nc, (0.0 OP s)); \ 194 MArray2<T> r (nr, nc, (0.0 OP s)); \
195 \ 195 \
196 for (int j = 0; j < nc; j++) \ 196 for (octave_idx_type j = 0; j < nc; j++) \
197 for (int i = a.cidx(j); i < a.cidx(j+1); i++) \ 197 for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++) \
198 r.elem (a.ridx (i), j) = a.data (i) OP s; \ 198 r.elem (a.ridx (i), j) = a.data (i) OP s; \
199 return r; \ 199 return r; \
200 } 200 }
201 201
202 #define SPARSE_A2S_OP_2(OP) \ 202 #define SPARSE_A2S_OP_2(OP) \
203 template <class T> \ 203 template <class T> \
204 MSparse<T> \ 204 MSparse<T> \
205 operator OP (const MSparse<T>& a, const T& s) \ 205 operator OP (const MSparse<T>& a, const T& s) \
206 { \ 206 { \
207 int nr = a.rows (); \ 207 octave_idx_type nr = a.rows (); \
208 int nc = a.cols (); \ 208 octave_idx_type nc = a.cols (); \
209 int nz = a.nnz (); \ 209 octave_idx_type nz = a.nnz (); \
210 \ 210 \
211 MSparse<T> r (nr, nc, nz); \ 211 MSparse<T> r (nr, nc, nz); \
212 \ 212 \
213 for (int i = 0; i < nz; i++) \ 213 for (octave_idx_type i = 0; i < nz; i++) \
214 { \ 214 { \
215 r.data(i) = a.data(i) OP s; \ 215 r.data(i) = a.data(i) OP s; \
216 r.ridx(i) = a.ridx(i); \ 216 r.ridx(i) = a.ridx(i); \
217 } \ 217 } \
218 for (int i = 0; i < nc + 1; i++) \ 218 for (octave_idx_type i = 0; i < nc + 1; i++) \
219 r.cidx(i) = a.cidx(i); \ 219 r.cidx(i) = a.cidx(i); \
220 r.maybe_compress (true); \ 220 r.maybe_compress (true); \
221 return r; \ 221 return r; \
222 } 222 }
223 223
232 #define SPARSE_SA2_OP_1(OP) \ 232 #define SPARSE_SA2_OP_1(OP) \
233 template <class T> \ 233 template <class T> \
234 MArray2<T> \ 234 MArray2<T> \
235 operator OP (const T& s, const MSparse<T>& a) \ 235 operator OP (const T& s, const MSparse<T>& a) \
236 { \ 236 { \
237 int nr = a.rows (); \ 237 octave_idx_type nr = a.rows (); \
238 int nc = a.cols (); \ 238 octave_idx_type nc = a.cols (); \
239 \ 239 \
240 MArray2<T> r (nr, nc, (s OP 0.0)); \ 240 MArray2<T> r (nr, nc, (s OP 0.0)); \
241 \ 241 \
242 for (int j = 0; j < nc; j++) \ 242 for (octave_idx_type j = 0; j < nc; j++) \
243 for (int i = a.cidx(j); i < a.cidx(j+1); i++) \ 243 for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++) \
244 r.elem (a.ridx (i), j) = s OP a.data (i); \ 244 r.elem (a.ridx (i), j) = s OP a.data (i); \
245 return r; \ 245 return r; \
246 } 246 }
247 247
248 #define SPARSE_SA2_OP_2(OP) \ 248 #define SPARSE_SA2_OP_2(OP) \
249 template <class T> \ 249 template <class T> \
250 MSparse<T> \ 250 MSparse<T> \
251 operator OP (const T& s, const MSparse<T>& a) \ 251 operator OP (const T& s, const MSparse<T>& a) \
252 { \ 252 { \
253 int nr = a.rows (); \ 253 octave_idx_type nr = a.rows (); \
254 int nc = a.cols (); \ 254 octave_idx_type nc = a.cols (); \
255 int nz = a.nnz (); \ 255 octave_idx_type nz = a.nnz (); \
256 \ 256 \
257 MSparse<T> r (nr, nc, nz); \ 257 MSparse<T> r (nr, nc, nz); \
258 \ 258 \
259 for (int i = 0; i < nz; i++) \ 259 for (octave_idx_type i = 0; i < nz; i++) \
260 { \ 260 { \
261 r.data(i) = s OP a.data(i); \ 261 r.data(i) = s OP a.data(i); \
262 r.ridx(i) = a.ridx(i); \ 262 r.ridx(i) = a.ridx(i); \
263 } \ 263 } \
264 for (int i = 0; i < nc + 1; i++) \ 264 for (octave_idx_type i = 0; i < nc + 1; i++) \
265 r.cidx(i) = a.cidx(i); \ 265 r.cidx(i) = a.cidx(i); \
266 r.maybe_compress (true); \ 266 r.maybe_compress (true); \
267 return r; \ 267 return r; \
268 } 268 }
269 269
279 MSparse<T> \ 279 MSparse<T> \
280 operator OP (const MSparse<T>& a, const MSparse<T>& b) \ 280 operator OP (const MSparse<T>& a, const MSparse<T>& b) \
281 { \ 281 { \
282 MSparse<T> r; \ 282 MSparse<T> r; \
283 \ 283 \
284 int a_nr = a.rows (); \ 284 octave_idx_type a_nr = a.rows (); \
285 int a_nc = a.cols (); \ 285 octave_idx_type a_nc = a.cols (); \
286 \ 286 \
287 int b_nr = b.rows (); \ 287 octave_idx_type b_nr = b.rows (); \
288 int b_nc = b.cols (); \ 288 octave_idx_type b_nc = b.cols (); \
289 \ 289 \
290 if (a_nr != b_nr || a_nc != b_nc) \ 290 if (a_nr != b_nr || a_nc != b_nc) \
291 gripe_nonconformant ("operator " # OP, a_nr, a_nc, b_nr, b_nc); \ 291 gripe_nonconformant ("operator " # OP, a_nr, a_nc, b_nr, b_nc); \
292 else \ 292 else \
293 { \ 293 { \
294 r = MSparse<T> (a_nr, a_nc, (a.nnz () + b.nnz ())); \ 294 r = MSparse<T> (a_nr, a_nc, (a.nnz () + b.nnz ())); \
295 \ 295 \
296 int jx = 0; \ 296 octave_idx_type jx = 0; \
297 r.cidx (0) = 0; \ 297 r.cidx (0) = 0; \
298 for (int i = 0 ; i < a_nc ; i++) \ 298 for (octave_idx_type i = 0 ; i < a_nc ; i++) \
299 { \ 299 { \
300 int ja = a.cidx(i); \ 300 octave_idx_type ja = a.cidx(i); \
301 int ja_max = a.cidx(i+1); \ 301 octave_idx_type ja_max = a.cidx(i+1); \
302 bool ja_lt_max= ja < ja_max; \ 302 bool ja_lt_max= ja < ja_max; \
303 \ 303 \
304 int jb = b.cidx(i); \ 304 octave_idx_type jb = b.cidx(i); \
305 int jb_max = b.cidx(i+1); \ 305 octave_idx_type jb_max = b.cidx(i+1); \
306 bool jb_lt_max = jb < jb_max; \ 306 bool jb_lt_max = jb < jb_max; \
307 \ 307 \
308 while (ja_lt_max || jb_lt_max ) \ 308 while (ja_lt_max || jb_lt_max ) \
309 { \ 309 { \
310 OCTAVE_QUIT; \ 310 OCTAVE_QUIT; \
354 MSparse<T> \ 354 MSparse<T> \
355 FCN (const MSparse<T>& a, const MSparse<T>& b) \ 355 FCN (const MSparse<T>& a, const MSparse<T>& b) \
356 { \ 356 { \
357 MSparse<T> r; \ 357 MSparse<T> r; \
358 \ 358 \
359 int a_nr = a.rows (); \ 359 octave_idx_type a_nr = a.rows (); \
360 int a_nc = a.cols (); \ 360 octave_idx_type a_nc = a.cols (); \
361 \ 361 \
362 int b_nr = b.rows (); \ 362 octave_idx_type b_nr = b.rows (); \
363 int b_nc = b.cols (); \ 363 octave_idx_type b_nc = b.cols (); \
364 \ 364 \
365 if (a_nr != b_nr || a_nc != b_nc) \ 365 if (a_nr != b_nr || a_nc != b_nc) \
366 gripe_nonconformant (#FCN, a_nr, a_nc, b_nr, b_nc); \ 366 gripe_nonconformant (#FCN, a_nr, a_nc, b_nr, b_nc); \
367 else \ 367 else \
368 { \ 368 { \
369 r = MSparse<T> (a_nr, a_nc, (a.nnz() > b.nnz() ? a.nnz() : b.nnz())); \ 369 r = MSparse<T> (a_nr, a_nc, (a.nnz() > b.nnz() ? a.nnz() : b.nnz())); \
370 \ 370 \
371 int jx = 0; \ 371 octave_idx_type jx = 0; \
372 r.cidx (0) = 0; \ 372 r.cidx (0) = 0; \
373 for (int i = 0 ; i < a_nc ; i++) \ 373 for (octave_idx_type i = 0 ; i < a_nc ; i++) \
374 { \ 374 { \
375 int ja = a.cidx(i); \ 375 octave_idx_type ja = a.cidx(i); \
376 int ja_max = a.cidx(i+1); \ 376 octave_idx_type ja_max = a.cidx(i+1); \
377 bool ja_lt_max= ja < ja_max; \ 377 bool ja_lt_max= ja < ja_max; \
378 \ 378 \
379 int jb = b.cidx(i); \ 379 octave_idx_type jb = b.cidx(i); \
380 int jb_max = b.cidx(i+1); \ 380 octave_idx_type jb_max = b.cidx(i+1); \
381 bool jb_lt_max = jb < jb_max; \ 381 bool jb_lt_max = jb < jb_max; \
382 \ 382 \
383 while (ja_lt_max || jb_lt_max ) \ 383 while (ja_lt_max || jb_lt_max ) \
384 { \ 384 { \
385 OCTAVE_QUIT; \ 385 OCTAVE_QUIT; \
420 FCN (const MSparse<T>& a, const MSparse<T>& b) \ 420 FCN (const MSparse<T>& a, const MSparse<T>& b) \
421 { \ 421 { \
422 MSparse<T> r; \ 422 MSparse<T> r; \
423 T Zero = T (); \ 423 T Zero = T (); \
424 \ 424 \
425 int a_nr = a.rows (); \ 425 octave_idx_type a_nr = a.rows (); \
426 int a_nc = a.cols (); \ 426 octave_idx_type a_nc = a.cols (); \
427 \ 427 \
428 int b_nr = b.rows (); \ 428 octave_idx_type b_nr = b.rows (); \
429 int b_nc = b.cols (); \ 429 octave_idx_type b_nc = b.cols (); \
430 \ 430 \
431 if (a_nr != b_nr || a_nc != b_nc) \ 431 if (a_nr != b_nr || a_nc != b_nc) \
432 gripe_nonconformant (#FCN, a_nr, a_nc, b_nr, b_nc); \ 432 gripe_nonconformant (#FCN, a_nr, a_nc, b_nr, b_nc); \
433 else \ 433 else \
434 { \ 434 { \
435 r = MSparse<T>( a_nr, a_nc, (Zero OP Zero)); \ 435 r = MSparse<T>( a_nr, a_nc, (Zero OP Zero)); \
436 \ 436 \
437 for (int i = 0 ; i < a_nc ; i++) \ 437 for (octave_idx_type i = 0 ; i < a_nc ; i++) \
438 { \ 438 { \
439 int ja = a.cidx(i); \ 439 octave_idx_type ja = a.cidx(i); \
440 int ja_max = a.cidx(i+1); \ 440 octave_idx_type ja_max = a.cidx(i+1); \
441 bool ja_lt_max= ja < ja_max; \ 441 bool ja_lt_max= ja < ja_max; \
442 \ 442 \
443 int jb = b.cidx(i); \ 443 octave_idx_type jb = b.cidx(i); \
444 int jb_max = b.cidx(i+1); \ 444 octave_idx_type jb_max = b.cidx(i+1); \
445 bool jb_lt_max = jb < jb_max; \ 445 bool jb_lt_max = jb < jb_max; \
446 \ 446 \
447 while (ja_lt_max || jb_lt_max ) \ 447 while (ja_lt_max || jb_lt_max ) \
448 { \ 448 { \
449 OCTAVE_QUIT; \ 449 OCTAVE_QUIT; \
491 template <class T> 491 template <class T>
492 MSparse<T> 492 MSparse<T>
493 operator - (const MSparse<T>& a) 493 operator - (const MSparse<T>& a)
494 { 494 {
495 MSparse<T> retval (a); 495 MSparse<T> retval (a);
496 int nz = a.nnz (); 496 octave_idx_type nz = a.nnz ();
497 for (int i = 0; i < nz; i++) 497 for (octave_idx_type i = 0; i < nz; i++)
498 retval.data(i) = - retval.data(i); 498 retval.data(i) = - retval.data(i);
499 return retval; 499 return retval;
500 } 500 }
501 501
502 /* 502 /*