comparison src/sparse-xpow.cc @ 5275:23b37da9fd5b

[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
55 octave_value 55 octave_value
56 xpow (const SparseMatrix& a, double b) 56 xpow (const SparseMatrix& a, double b)
57 { 57 {
58 octave_value retval; 58 octave_value retval;
59 59
60 int nr = a.rows (); 60 octave_idx_type nr = a.rows ();
61 int nc = a.cols (); 61 octave_idx_type nc = a.cols ();
62 62
63 if (nr == 0 || nc == 0 || nr != nc) 63 if (nr == 0 || nc == 0 || nr != nc)
64 error ("for A^b, A must be square"); 64 error ("for A^b, A must be square");
65 else 65 else
66 { 66 {
68 { 68 {
69 int btmp = static_cast<int> (b); 69 int btmp = static_cast<int> (b);
70 if (btmp == 0) 70 if (btmp == 0)
71 { 71 {
72 SparseMatrix tmp = SparseMatrix (nr, nr, nr); 72 SparseMatrix tmp = SparseMatrix (nr, nr, nr);
73 for (int i = 0; i < nr; i++) 73 for (octave_idx_type i = 0; i < nr; i++)
74 { 74 {
75 tmp.data (i) = 1.0; 75 tmp.data (i) = 1.0;
76 tmp.ridx (i) = i; 76 tmp.ridx (i) = i;
77 } 77 }
78 for (int i = 0; i < nr + 1; i++) 78 for (octave_idx_type i = 0; i < nr + 1; i++)
79 tmp.cidx (i) = i; 79 tmp.cidx (i) = i;
80 80
81 retval = tmp; 81 retval = tmp;
82 } 82 }
83 else 83 else
85 SparseMatrix atmp; 85 SparseMatrix atmp;
86 if (btmp < 0) 86 if (btmp < 0)
87 { 87 {
88 btmp = -btmp; 88 btmp = -btmp;
89 89
90 int info; 90 octave_idx_type info;
91 double rcond = 0.0; 91 double rcond = 0.0;
92 92
93 atmp = a.inverse (info, rcond, 1); 93 atmp = a.inverse (info, rcond, 1);
94 94
95 if (info == -1) 95 if (info == -1)
127 octave_value 127 octave_value
128 xpow (const SparseComplexMatrix& a, double b) 128 xpow (const SparseComplexMatrix& a, double b)
129 { 129 {
130 octave_value retval; 130 octave_value retval;
131 131
132 int nr = a.rows (); 132 octave_idx_type nr = a.rows ();
133 int nc = a.cols (); 133 octave_idx_type nc = a.cols ();
134 134
135 if (nr == 0 || nc == 0 || nr != nc) 135 if (nr == 0 || nc == 0 || nr != nc)
136 error ("for A^b, A must be square"); 136 error ("for A^b, A must be square");
137 else 137 else
138 { 138 {
140 { 140 {
141 int btmp = static_cast<int> (b); 141 int btmp = static_cast<int> (b);
142 if (btmp == 0) 142 if (btmp == 0)
143 { 143 {
144 SparseMatrix tmp = SparseMatrix (nr, nr, nr); 144 SparseMatrix tmp = SparseMatrix (nr, nr, nr);
145 for (int i = 0; i < nr; i++) 145 for (octave_idx_type i = 0; i < nr; i++)
146 { 146 {
147 tmp.data (i) = 1.0; 147 tmp.data (i) = 1.0;
148 tmp.ridx (i) = i; 148 tmp.ridx (i) = i;
149 } 149 }
150 for (int i = 0; i < nr + 1; i++) 150 for (octave_idx_type i = 0; i < nr + 1; i++)
151 tmp.cidx (i) = i; 151 tmp.cidx (i) = i;
152 152
153 retval = tmp; 153 retval = tmp;
154 } 154 }
155 else 155 else
157 SparseComplexMatrix atmp; 157 SparseComplexMatrix atmp;
158 if (btmp < 0) 158 if (btmp < 0)
159 { 159 {
160 btmp = -btmp; 160 btmp = -btmp;
161 161
162 int info; 162 octave_idx_type info;
163 double rcond = 0.0; 163 double rcond = 0.0;
164 164
165 atmp = a.inverse (info, rcond, 1); 165 atmp = a.inverse (info, rcond, 1);
166 166
167 if (info == -1) 167 if (info == -1)
229 octave_value 229 octave_value
230 elem_xpow (double a, const SparseMatrix& b) 230 elem_xpow (double a, const SparseMatrix& b)
231 { 231 {
232 octave_value retval; 232 octave_value retval;
233 233
234 int nr = b.rows (); 234 octave_idx_type nr = b.rows ();
235 int nc = b.cols (); 235 octave_idx_type nc = b.cols ();
236 236
237 double d1, d2; 237 double d1, d2;
238 238
239 if (a < 0.0 && ! b.all_integers (d1, d2)) 239 if (a < 0.0 && ! b.all_integers (d1, d2))
240 { 240 {
241 Complex atmp (a); 241 Complex atmp (a);
242 ComplexMatrix result (nr, nc); 242 ComplexMatrix result (nr, nc);
243 243
244 for (int j = 0; j < nc; j++) 244 for (octave_idx_type j = 0; j < nc; j++)
245 { 245 {
246 for (int i = 0; i < nr; i++) 246 for (octave_idx_type i = 0; i < nr; i++)
247 { 247 {
248 OCTAVE_QUIT; 248 OCTAVE_QUIT;
249 result (i, j) = pow (atmp, b(i,j)); 249 result (i, j) = pow (atmp, b(i,j));
250 } 250 }
251 } 251 }
254 } 254 }
255 else 255 else
256 { 256 {
257 Matrix result (nr, nc); 257 Matrix result (nr, nc);
258 258
259 for (int j = 0; j < nc; j++) 259 for (octave_idx_type j = 0; j < nc; j++)
260 { 260 {
261 for (int i = 0; i < nr; i++) 261 for (octave_idx_type i = 0; i < nr; i++)
262 { 262 {
263 OCTAVE_QUIT; 263 OCTAVE_QUIT;
264 result (i, j) = pow (a, b(i,j)); 264 result (i, j) = pow (a, b(i,j));
265 } 265 }
266 } 266 }
273 273
274 // -*- 2 -*- 274 // -*- 2 -*-
275 octave_value 275 octave_value
276 elem_xpow (double a, const SparseComplexMatrix& b) 276 elem_xpow (double a, const SparseComplexMatrix& b)
277 { 277 {
278 int nr = b.rows (); 278 octave_idx_type nr = b.rows ();
279 int nc = b.cols (); 279 octave_idx_type nc = b.cols ();
280 280
281 Complex atmp (a); 281 Complex atmp (a);
282 ComplexMatrix result (nr, nc); 282 ComplexMatrix result (nr, nc);
283 283
284 for (int j = 0; j < nc; j++) 284 for (octave_idx_type j = 0; j < nc; j++)
285 { 285 {
286 for (int i = 0; i < nr; i++) 286 for (octave_idx_type i = 0; i < nr; i++)
287 { 287 {
288 OCTAVE_QUIT; 288 OCTAVE_QUIT;
289 result (i, j) = pow (atmp, b(i,j)); 289 result (i, j) = pow (atmp, b(i,j));
290 } 290 }
291 } 291 }
301 // sparse matrix with same structure as a, which is strictly 301 // sparse matrix with same structure as a, which is strictly
302 // incorrect. Keep compatiability. 302 // incorrect. Keep compatiability.
303 303
304 octave_value retval; 304 octave_value retval;
305 305
306 int nz = a.nnz (); 306 octave_idx_type nz = a.nnz ();
307 307
308 if (b <= 0.0) 308 if (b <= 0.0)
309 { 309 {
310 int nr = a.rows (); 310 octave_idx_type nr = a.rows ();
311 int nc = a.cols (); 311 octave_idx_type nc = a.cols ();
312 312
313 if (static_cast<int> (b) != b && a.any_element_is_negative ()) 313 if (static_cast<int> (b) != b && a.any_element_is_negative ())
314 { 314 {
315 ComplexMatrix result (nr, nc, Complex (pow (0.0, b))); 315 ComplexMatrix result (nr, nc, Complex (pow (0.0, b)));
316 316
317 // XXX FIXME XXX -- avoid apparent GNU libm bug by 317 // XXX FIXME XXX -- avoid apparent GNU libm bug by
318 // converting A and B to complex instead of just A. 318 // converting A and B to complex instead of just A.
319 Complex btmp (b); 319 Complex btmp (b);
320 320
321 for (int j = 0; j < nc; j++) 321 for (octave_idx_type j = 0; j < nc; j++)
322 for (int i = a.cidx(j); i < a.cidx(j+1); i++) 322 for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
323 { 323 {
324 OCTAVE_QUIT; 324 OCTAVE_QUIT;
325 325
326 Complex atmp (a.data (i)); 326 Complex atmp (a.data (i));
327 327
332 } 332 }
333 else 333 else
334 { 334 {
335 Matrix result (nr, nc, (pow (0.0, b))); 335 Matrix result (nr, nc, (pow (0.0, b)));
336 336
337 for (int j = 0; j < nc; j++) 337 for (octave_idx_type j = 0; j < nc; j++)
338 for (int i = a.cidx(j); i < a.cidx(j+1); i++) 338 for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
339 { 339 {
340 OCTAVE_QUIT; 340 OCTAVE_QUIT;
341 result (a.ridx(i), j) = pow (a.data (i), b); 341 result (a.ridx(i), j) = pow (a.data (i), b);
342 } 342 }
343 343
346 } 346 }
347 else if (static_cast<int> (b) != b && a.any_element_is_negative ()) 347 else if (static_cast<int> (b) != b && a.any_element_is_negative ())
348 { 348 {
349 SparseComplexMatrix result (a); 349 SparseComplexMatrix result (a);
350 350
351 for (int i = 0; i < nz; i++) 351 for (octave_idx_type i = 0; i < nz; i++)
352 { 352 {
353 OCTAVE_QUIT; 353 OCTAVE_QUIT;
354 354
355 // XXX FIXME XXX -- avoid apparent GNU libm bug by 355 // XXX FIXME XXX -- avoid apparent GNU libm bug by
356 // converting A and B to complex instead of just A. 356 // converting A and B to complex instead of just A.
367 } 367 }
368 else 368 else
369 { 369 {
370 SparseMatrix result (a); 370 SparseMatrix result (a);
371 371
372 for (int i = 0; i < nz; i++) 372 for (octave_idx_type i = 0; i < nz; i++)
373 { 373 {
374 OCTAVE_QUIT; 374 OCTAVE_QUIT;
375 result.data (i) = pow (a.data (i), b); 375 result.data (i) = pow (a.data (i), b);
376 } 376 }
377 377
387 octave_value 387 octave_value
388 elem_xpow (const SparseMatrix& a, const SparseMatrix& b) 388 elem_xpow (const SparseMatrix& a, const SparseMatrix& b)
389 { 389 {
390 octave_value retval; 390 octave_value retval;
391 391
392 int nr = a.rows (); 392 octave_idx_type nr = a.rows ();
393 int nc = a.cols (); 393 octave_idx_type nc = a.cols ();
394 394
395 int b_nr = b.rows (); 395 octave_idx_type b_nr = b.rows ();
396 int b_nc = b.cols (); 396 octave_idx_type b_nc = b.cols ();
397 397
398 if (nr != b_nr || nc != b_nc) 398 if (nr != b_nr || nc != b_nc)
399 { 399 {
400 gripe_nonconformant ("operator .^", nr, nc, b_nr, b_nc); 400 gripe_nonconformant ("operator .^", nr, nc, b_nr, b_nc);
401 return octave_value (); 401 return octave_value ();
402 } 402 }
403 403
404 int convert_to_complex = 0; 404 int convert_to_complex = 0;
405 for (int j = 0; j < nc; j++) 405 for (octave_idx_type j = 0; j < nc; j++)
406 for (int i = 0; i < nr; i++) 406 for (octave_idx_type i = 0; i < nr; i++)
407 { 407 {
408 OCTAVE_QUIT; 408 OCTAVE_QUIT;
409 double atmp = a (i, j); 409 double atmp = a (i, j);
410 double btmp = b (i, j); 410 double btmp = b (i, j);
411 if (atmp < 0.0 && static_cast<int> (btmp) != btmp) 411 if (atmp < 0.0 && static_cast<int> (btmp) != btmp)
415 } 415 }
416 } 416 }
417 417
418 done: 418 done:
419 419
420 int nel = 0; 420 octave_idx_type nel = 0;
421 for (int j = 0; j < nc; j++) 421 for (octave_idx_type j = 0; j < nc; j++)
422 for (int i = 0; i < nr; i++) 422 for (octave_idx_type i = 0; i < nr; i++)
423 if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.)) 423 if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.))
424 nel++; 424 nel++;
425 425
426 if (convert_to_complex) 426 if (convert_to_complex)
427 { 427 {
428 SparseComplexMatrix complex_result (nr, nc, nel); 428 SparseComplexMatrix complex_result (nr, nc, nel);
429 429
430 int ii = 0; 430 octave_idx_type ii = 0;
431 complex_result.cidx(0) = 0; 431 complex_result.cidx(0) = 0;
432 for (int j = 0; j < nc; j++) 432 for (octave_idx_type j = 0; j < nc; j++)
433 { 433 {
434 for (int i = 0; i < nr; i++) 434 for (octave_idx_type i = 0; i < nr; i++)
435 { 435 {
436 OCTAVE_QUIT; 436 OCTAVE_QUIT;
437 Complex atmp (a (i, j)); 437 Complex atmp (a (i, j));
438 Complex btmp (b (i, j)); 438 Complex btmp (b (i, j));
439 Complex tmp = pow (atmp, btmp); 439 Complex tmp = pow (atmp, btmp);
450 retval = complex_result; 450 retval = complex_result;
451 } 451 }
452 else 452 else
453 { 453 {
454 SparseMatrix result (nr, nc, nel); 454 SparseMatrix result (nr, nc, nel);
455 int ii = 0; 455 octave_idx_type ii = 0;
456 456
457 result.cidx (0) = 0; 457 result.cidx (0) = 0;
458 for (int j = 0; j < nc; j++) 458 for (octave_idx_type j = 0; j < nc; j++)
459 { 459 {
460 for (int i = 0; i < nr; i++) 460 for (octave_idx_type i = 0; i < nr; i++)
461 { 461 {
462 OCTAVE_QUIT; 462 OCTAVE_QUIT;
463 double tmp = pow (a (i, j), b (i, j)); 463 double tmp = pow (a (i, j), b (i, j));
464 if (tmp != 0.) 464 if (tmp != 0.)
465 { 465 {
487 if (b == 0.0) 487 if (b == 0.0)
488 // Can this case ever happen, due to automatic retyping with maybe_mutate? 488 // Can this case ever happen, due to automatic retyping with maybe_mutate?
489 retval = octave_value (NDArray (a.dims (), 1)); 489 retval = octave_value (NDArray (a.dims (), 1));
490 else 490 else
491 { 491 {
492 int nz = a.nnz (); 492 octave_idx_type nz = a.nnz ();
493 SparseComplexMatrix result (a); 493 SparseComplexMatrix result (a);
494 494
495 for (int i = 0; i < nz; i++) 495 for (octave_idx_type i = 0; i < nz; i++)
496 { 496 {
497 OCTAVE_QUIT; 497 OCTAVE_QUIT;
498 result.data (i) = pow (Complex (a.data (i)), b); 498 result.data (i) = pow (Complex (a.data (i)), b);
499 } 499 }
500 500
508 508
509 // -*- 6 -*- 509 // -*- 6 -*-
510 octave_value 510 octave_value
511 elem_xpow (const SparseMatrix& a, const SparseComplexMatrix& b) 511 elem_xpow (const SparseMatrix& a, const SparseComplexMatrix& b)
512 { 512 {
513 int nr = a.rows (); 513 octave_idx_type nr = a.rows ();
514 int nc = a.cols (); 514 octave_idx_type nc = a.cols ();
515 515
516 int b_nr = b.rows (); 516 octave_idx_type b_nr = b.rows ();
517 int b_nc = b.cols (); 517 octave_idx_type b_nc = b.cols ();
518 518
519 if (nr != b_nr || nc != b_nc) 519 if (nr != b_nr || nc != b_nc)
520 { 520 {
521 gripe_nonconformant ("operator .^", nr, nc, b_nr, b_nc); 521 gripe_nonconformant ("operator .^", nr, nc, b_nr, b_nc);
522 return octave_value (); 522 return octave_value ();
523 } 523 }
524 524
525 int nel = 0; 525 octave_idx_type nel = 0;
526 for (int j = 0; j < nc; j++) 526 for (octave_idx_type j = 0; j < nc; j++)
527 for (int i = 0; i < nr; i++) 527 for (octave_idx_type i = 0; i < nr; i++)
528 if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.)) 528 if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.))
529 nel++; 529 nel++;
530 530
531 SparseComplexMatrix result (nr, nc, nel); 531 SparseComplexMatrix result (nr, nc, nel);
532 int ii = 0; 532 octave_idx_type ii = 0;
533 533
534 result.cidx(0) = 0; 534 result.cidx(0) = 0;
535 for (int j = 0; j < nc; j++) 535 for (octave_idx_type j = 0; j < nc; j++)
536 { 536 {
537 for (int i = 0; i < nr; i++) 537 for (octave_idx_type i = 0; i < nr; i++)
538 { 538 {
539 OCTAVE_QUIT; 539 OCTAVE_QUIT;
540 Complex tmp = pow (Complex (a (i, j)), b (i, j)); 540 Complex tmp = pow (Complex (a (i, j)), b (i, j));
541 if (tmp != 0.) 541 if (tmp != 0.)
542 { 542 {
554 554
555 // -*- 7 -*- 555 // -*- 7 -*-
556 octave_value 556 octave_value
557 elem_xpow (const Complex& a, const SparseMatrix& b) 557 elem_xpow (const Complex& a, const SparseMatrix& b)
558 { 558 {
559 int nr = b.rows (); 559 octave_idx_type nr = b.rows ();
560 int nc = b.cols (); 560 octave_idx_type nc = b.cols ();
561 561
562 ComplexMatrix result (nr, nc); 562 ComplexMatrix result (nr, nc);
563 563
564 for (int j = 0; j < nc; j++) 564 for (octave_idx_type j = 0; j < nc; j++)
565 { 565 {
566 for (int i = 0; i < nr; i++) 566 for (octave_idx_type i = 0; i < nr; i++)
567 { 567 {
568 OCTAVE_QUIT; 568 OCTAVE_QUIT;
569 double btmp = b (i, j); 569 double btmp = b (i, j);
570 if (xisint (btmp)) 570 if (xisint (btmp))
571 result (i, j) = pow (a, static_cast<int> (btmp)); 571 result (i, j) = pow (a, static_cast<int> (btmp));
579 579
580 // -*- 8 -*- 580 // -*- 8 -*-
581 octave_value 581 octave_value
582 elem_xpow (const Complex& a, const SparseComplexMatrix& b) 582 elem_xpow (const Complex& a, const SparseComplexMatrix& b)
583 { 583 {
584 int nr = b.rows (); 584 octave_idx_type nr = b.rows ();
585 int nc = b.cols (); 585 octave_idx_type nc = b.cols ();
586 586
587 ComplexMatrix result (nr, nc); 587 ComplexMatrix result (nr, nc);
588 for (int j = 0; j < nc; j++) 588 for (octave_idx_type j = 0; j < nc; j++)
589 for (int i = 0; i < nr; i++) 589 for (octave_idx_type i = 0; i < nr; i++)
590 { 590 {
591 OCTAVE_QUIT; 591 OCTAVE_QUIT;
592 result (i, j) = pow (a, b (i, j)); 592 result (i, j) = pow (a, b (i, j));
593 } 593 }
594 594
601 { 601 {
602 octave_value retval; 602 octave_value retval;
603 603
604 if (b <= 0) 604 if (b <= 0)
605 { 605 {
606 int nr = a.rows (); 606 octave_idx_type nr = a.rows ();
607 int nc = a.cols (); 607 octave_idx_type nc = a.cols ();
608 608
609 ComplexMatrix result (nr, nc, Complex (pow (0.0, b))); 609 ComplexMatrix result (nr, nc, Complex (pow (0.0, b)));
610 610
611 if (xisint (b)) 611 if (xisint (b))
612 { 612 {
613 for (int j = 0; j < nc; j++) 613 for (octave_idx_type j = 0; j < nc; j++)
614 for (int i = a.cidx(j); i < a.cidx(j+1); i++) 614 for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
615 { 615 {
616 OCTAVE_QUIT; 616 OCTAVE_QUIT;
617 result (a.ridx(i), j) = 617 result (a.ridx(i), j) =
618 pow (a.data (i), static_cast<int> (b)); 618 pow (a.data (i), static_cast<int> (b));
619 } 619 }
620 } 620 }
621 else 621 else
622 { 622 {
623 for (int j = 0; j < nc; j++) 623 for (octave_idx_type j = 0; j < nc; j++)
624 for (int i = a.cidx(j); i < a.cidx(j+1); i++) 624 for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
625 { 625 {
626 OCTAVE_QUIT; 626 OCTAVE_QUIT;
627 result (a.ridx(i), j) = pow (a.data (i), b); 627 result (a.ridx(i), j) = pow (a.data (i), b);
628 } 628 }
629 } 629 }
630 630
631 retval = result; 631 retval = result;
632 } 632 }
633 else 633 else
634 { 634 {
635 int nz = a.nnz (); 635 octave_idx_type nz = a.nnz ();
636 636
637 SparseComplexMatrix result (a); 637 SparseComplexMatrix result (a);
638 638
639 if (xisint (b)) 639 if (xisint (b))
640 { 640 {
641 for (int i = 0; i < nz; i++) 641 for (octave_idx_type i = 0; i < nz; i++)
642 { 642 {
643 OCTAVE_QUIT; 643 OCTAVE_QUIT;
644 result.data (i) = pow (a.data (i), static_cast<int> (b)); 644 result.data (i) = pow (a.data (i), static_cast<int> (b));
645 } 645 }
646 } 646 }
647 else 647 else
648 { 648 {
649 for (int i = 0; i < nz; i++) 649 for (octave_idx_type i = 0; i < nz; i++)
650 { 650 {
651 OCTAVE_QUIT; 651 OCTAVE_QUIT;
652 result.data (i) = pow (a.data (i), b); 652 result.data (i) = pow (a.data (i), b);
653 } 653 }
654 } 654 }
663 663
664 // -*- 10 -*- 664 // -*- 10 -*-
665 octave_value 665 octave_value
666 elem_xpow (const SparseComplexMatrix& a, const SparseMatrix& b) 666 elem_xpow (const SparseComplexMatrix& a, const SparseMatrix& b)
667 { 667 {
668 int nr = a.rows (); 668 octave_idx_type nr = a.rows ();
669 int nc = a.cols (); 669 octave_idx_type nc = a.cols ();
670 670
671 int b_nr = b.rows (); 671 octave_idx_type b_nr = b.rows ();
672 int b_nc = b.cols (); 672 octave_idx_type b_nc = b.cols ();
673 673
674 if (nr != b_nr || nc != b_nc) 674 if (nr != b_nr || nc != b_nc)
675 { 675 {
676 gripe_nonconformant ("operator .^", nr, nc, b_nr, b_nc); 676 gripe_nonconformant ("operator .^", nr, nc, b_nr, b_nc);
677 return octave_value (); 677 return octave_value ();
678 } 678 }
679 679
680 int nel = 0; 680 octave_idx_type nel = 0;
681 for (int j = 0; j < nc; j++) 681 for (octave_idx_type j = 0; j < nc; j++)
682 for (int i = 0; i < nr; i++) 682 for (octave_idx_type i = 0; i < nr; i++)
683 if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.)) 683 if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.))
684 nel++; 684 nel++;
685 685
686 SparseComplexMatrix result (nr, nc, nel); 686 SparseComplexMatrix result (nr, nc, nel);
687 int ii = 0; 687 octave_idx_type ii = 0;
688 688
689 result.cidx (0) = 0; 689 result.cidx (0) = 0;
690 for (int j = 0; j < nc; j++) 690 for (octave_idx_type j = 0; j < nc; j++)
691 { 691 {
692 for (int i = 0; i < nr; i++) 692 for (octave_idx_type i = 0; i < nr; i++)
693 { 693 {
694 OCTAVE_QUIT; 694 OCTAVE_QUIT;
695 double btmp = b (i, j); 695 double btmp = b (i, j);
696 Complex tmp; 696 Complex tmp;
697 697
723 // Can this case ever happen, due to automatic retyping with maybe_mutate? 723 // Can this case ever happen, due to automatic retyping with maybe_mutate?
724 retval = octave_value (NDArray (a.dims (), 1)); 724 retval = octave_value (NDArray (a.dims (), 1));
725 else 725 else
726 { 726 {
727 727
728 int nz = a.nnz (); 728 octave_idx_type nz = a.nnz ();
729 729
730 SparseComplexMatrix result (a); 730 SparseComplexMatrix result (a);
731 731
732 for (int i = 0; i < nz; i++) 732 for (octave_idx_type i = 0; i < nz; i++)
733 { 733 {
734 OCTAVE_QUIT; 734 OCTAVE_QUIT;
735 result.data (i) = pow (a.data (i), b); 735 result.data (i) = pow (a.data (i), b);
736 } 736 }
737 737
745 745
746 // -*- 12 -*- 746 // -*- 12 -*-
747 octave_value 747 octave_value
748 elem_xpow (const SparseComplexMatrix& a, const SparseComplexMatrix& b) 748 elem_xpow (const SparseComplexMatrix& a, const SparseComplexMatrix& b)
749 { 749 {
750 int nr = a.rows (); 750 octave_idx_type nr = a.rows ();
751 int nc = a.cols (); 751 octave_idx_type nc = a.cols ();
752 752
753 int b_nr = b.rows (); 753 octave_idx_type b_nr = b.rows ();
754 int b_nc = b.cols (); 754 octave_idx_type b_nc = b.cols ();
755 755
756 if (nr != b_nr || nc != b_nc) 756 if (nr != b_nr || nc != b_nc)
757 { 757 {
758 gripe_nonconformant ("operator .^", nr, nc, b_nr, b_nc); 758 gripe_nonconformant ("operator .^", nr, nc, b_nr, b_nc);
759 return octave_value (); 759 return octave_value ();
760 } 760 }
761 761
762 int nel = 0; 762 octave_idx_type nel = 0;
763 for (int j = 0; j < nc; j++) 763 for (octave_idx_type j = 0; j < nc; j++)
764 for (int i = 0; i < nr; i++) 764 for (octave_idx_type i = 0; i < nr; i++)
765 if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.)) 765 if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.))
766 nel++; 766 nel++;
767 767
768 SparseComplexMatrix result (nr, nc, nel); 768 SparseComplexMatrix result (nr, nc, nel);
769 int ii = 0; 769 octave_idx_type ii = 0;
770 770
771 result.cidx (0) = 0; 771 result.cidx (0) = 0;
772 for (int j = 0; j < nc; j++) 772 for (octave_idx_type j = 0; j < nc; j++)
773 { 773 {
774 for (int i = 0; i < nr; i++) 774 for (octave_idx_type i = 0; i < nr; i++)
775 { 775 {
776 OCTAVE_QUIT; 776 OCTAVE_QUIT;
777 Complex tmp = pow (a (i, j), b (i, j)); 777 Complex tmp = pow (a (i, j), b (i, j));
778 if (tmp != 0.) 778 if (tmp != 0.)
779 { 779 {