C-XSC - A C++ Class Library for Extended Scientific Computing 2.5.4
scmatrix.hpp
1/*
2** CXSC is a C++ library for eXtended Scientific Computing (V 2.5.4)
3**
4** Copyright (C) 1990-2000 Institut fuer Angewandte Mathematik,
5** Universitaet Karlsruhe, Germany
6** (C) 2000-2014 Wiss. Rechnen/Softwaretechnologie
7** Universitaet Wuppertal, Germany
8**
9** This library is free software; you can redistribute it and/or
10** modify it under the terms of the GNU Library General Public
11** License as published by the Free Software Foundation; either
12** version 2 of the License, or (at your option) any later version.
13**
14** This library is distributed in the hope that it will be useful,
15** but WITHOUT ANY WARRANTY; without even the implied warranty of
16** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17** Library General Public License for more details.
18**
19** You should have received a copy of the GNU Library General Public
20** License along with this library; if not, write to the Free
21** Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
22*/
23
24/* CVS $Id: scmatrix.hpp,v 1.20 2014/01/30 17:23:48 cxsc Exp $ */
25
26#ifndef _CXSC_SCMATRIX_HPP_INCLUDED
27#define _CXSC_SCMATRIX_HPP_INCLUDED
28
29#include <complex.hpp>
30#include <cmatrix.hpp>
31#include <scvector.hpp>
32#include <cidot.hpp>
33#include <vector>
34#include <algorithm>
35#include <iostream>
36#include <sparsecdot.hpp>
37#include <sparsematrix.hpp>
38#include <srmatrix.hpp>
39
40namespace cxsc {
41
42//definiert in srmatrix.hpp
43//enum STORAGE_TYPE{triplet,compressed_row,compressed_column};
44
45class scmatrix_slice;
46class scmatrix_subv;
47class scimatrix;
48class scimatrix_slice;
49class scimatrix_subv;
50
51
52inline bool comp_pair_c(std::pair<int,complex> p1, std::pair<int,complex> p2) {
53 return p1.first < p2.first;
54}
55
57
69class scmatrix {
70
71 private:
72 std::vector<int> p;
73 std::vector<int> ind;
74 std::vector<complex> x;
75 int m;
76 int n;
77 int lb1,ub1,lb2,ub2;
78
79 public:
80
82 std::vector<int>& column_pointers() {
83 return p;
84 }
85
87 std::vector<int>& row_indices() {
88 return ind;
89 }
90
92 std::vector<complex>& values() {
93 return x;
94 }
95
97 const std::vector<int>& column_pointers() const {
98 return p;
99 }
100
102 const std::vector<int>& row_indices() const {
103 return ind;
104 }
105
107 const std::vector<complex>& values() const {
108 return x;
109 }
110
113 p.push_back(0);
114 m = n = 0;
115 lb1 = lb2 = ub1 = ub2 = 0;
116 }
117
119 scmatrix(const int r, const int c) : m(r),n(c),lb1(1),ub1(r),lb2(1),ub2(c) {
120 p = std::vector<int>((n>0) ? n+1 : 1, 0);
121 ind.reserve(2*(m+n));
122 x.reserve(2*(m+n));
123
124 p[0] = 0;
125 }
126
128 scmatrix(const int r, const int c, const int e) : m(r),n(c),lb1(1),ub1(r),lb2(1),ub2(c) {
129 p = std::vector<int>((n>0) ? n+1 : 1, 0);
130 ind.reserve(e);
131 x.reserve(e);
132
133 p[0] = 0;
134 }
135
137
143 scmatrix(const int m, const int n, const int nnz, const intvector& rows, const intvector& cols, const cvector& values, const enum STORAGE_TYPE t = triplet) {
144 if(t == triplet) {
145 this->m = m;
146 this->n = n;
147 p = std::vector<int>(n+1,0);
148 ind.reserve(nnz);
149 x.reserve(nnz);
150 lb1 = lb2 = 1;
151 ub1 = m; ub2 = n;
152
153 std::vector<triplet_store<complex> > work;
154 work.reserve(nnz);
155
156 for(int k=0 ; k<nnz ; k++) {
157 work.push_back(triplet_store<complex>(rows[Lb(rows)+k],cols[Lb(cols)+k],values[Lb(values)+k]));
158 }
159
160 sort(work.begin(), work.end());
161
162 int i=0;
163
164 for(int j=0 ; j<n ; j++) {
165
166 while((unsigned int)i < work.size() && work[i].col == j ) {
167 ind.push_back(work[i].row);
168 x.push_back(work[i].val);
169 i++;
170 }
171
172 p[j+1] = i;
173 }
174
175 } else if(t == compressed_row) {
176
177 this->m = m;
178 this->n = n;
179 p = std::vector<int>(n+1,0);
180 ind.reserve(nnz);
181 x.reserve(nnz);
182 lb1 = lb2 = 1;
183 ub1 = m; ub2 = n;
184
185 for(int i=0 ; i<n+1 ; i++)
186 p[i] = rows[Lb(rows)+i];
187
188 std::vector<triplet_store<complex> > work;
189 work.reserve(nnz);
190
191 for(int j=0 ; j<n ; j++) {
192 for(int k=p[j] ; k<p[j+1] ; k++) {
193 work.push_back(triplet_store<complex>(j,cols[Lb(cols)+k],values[Lb(values)+k]));
194 }
195 }
196
197 sort(work.begin(), work.end());
198
199 int i=0;
200
201 for(int j=0 ; j<n ; j++) {
202
203 while((unsigned int)i < work.size() && work[i].col == j ) {
204 ind.push_back(work[i].row);
205 x.push_back(work[i].val);
206 i++;
207 }
208
209 p[j+1] = i;
210 }
211
212 } else if(t == compressed_column) {
213 this->m = m;
214 this->n = n;
215 p = std::vector<int>(n+1,0);
216 ind.reserve(nnz);
217 x.reserve(nnz);
218 lb1 = lb2 = 1;
219 ub1 = m; ub2 = n;
220
221 for(int i=0 ; i<n+1 ; i++)
222 p[i] = rows[Lb(rows)+i];
223
224 std::vector<std::pair<int,complex> > work;
225 work.reserve(n);
226
227 for(int j=0 ; j<n ; j++) {
228 work.clear();
229
230 for(int k=p[j] ; k<p[j+1] ; k++) {
231 work.push_back(std::make_pair(cols[Lb(cols)+k],values[Lb(values)+k]));
232 }
233
234 std::sort(work.begin(),work.end(),comp_pair_c);
235
236 for(unsigned int i=0 ; i<work.size() ; i++) {
237 ind.push_back(work[i].first);
238 x.push_back(work[i].second);
239 }
240 }
241
242 }
243
244 }
245
247
254 scmatrix(const int m, const int n, const int nnz, const int* rows, const int* cols, const complex* values, const enum STORAGE_TYPE t = triplet) {
255 if(t == triplet) {
256 this->m = m;
257 this->n = n;
258 p = std::vector<int>(n+1,0);
259 ind.reserve(nnz);
260 x.reserve(nnz);
261 lb1 = lb2 = 1;
262 ub1 = m; ub2 = n;
263
264 std::vector<triplet_store<complex> > work;
265 work.reserve(nnz);
266
267 for(int k=0 ; k<nnz ; k++) {
268 work.push_back(triplet_store<complex>(rows[k],cols[k],values[k]));
269 }
270
271 sort(work.begin(), work.end());
272
273 int i=0;
274
275 for(int j=0 ; j<n ; j++) {
276
277 while((unsigned int)i < work.size() && work[i].col == j ) {
278 ind.push_back(work[i].row);
279 x.push_back(work[i].val);
280 i++;
281 }
282
283 p[j+1] = i;
284 }
285
286 } else if(t == compressed_row) {
287
288 this->m = m;
289 this->n = n;
290 p = std::vector<int>(n+1,0);
291 ind.reserve(nnz);
292 x.reserve(nnz);
293 lb1 = lb2 = 1;
294 ub1 = m; ub2 = n;
295
296 for(int i=0 ; i<n+1 ; i++)
297 p[i] = rows[i];
298
299 std::vector<triplet_store<complex> > work;
300 work.reserve(nnz);
301
302 for(int j=0 ; j<n ; j++) {
303 for(int k=p[j] ; k<p[j+1] ; k++) {
304 work.push_back(triplet_store<complex>(j,cols[k],values[k]));
305 }
306 }
307
308 sort(work.begin(), work.end());
309
310 int i=0;
311
312 for(int j=0 ; j<n ; j++) {
313
314 while((unsigned int)i < work.size() && work[i].col == j ) {
315 ind.push_back(work[i].row);
316 x.push_back(work[i].val);
317 i++;
318 }
319
320 p[j+1] = i;
321 }
322
323 } else if(t == compressed_column) {
324 this->m = m;
325 this->n = n;
326 p = std::vector<int>(n+1,0);
327 ind.reserve(nnz);
328 x.reserve(nnz);
329 lb1 = lb2 = 1;
330 ub1 = m; ub2 = n;
331
332 for(int i=0 ; i<n+1 ; i++)
333 p[i] = rows[i];
334
335 std::vector<std::pair<int,complex> > work;
336 work.reserve(n);
337
338 for(int j=0 ; j<n ; j++) {
339 work.clear();
340
341 for(int k=p[j] ; k<p[j+1] ; k++) {
342 work.push_back(std::make_pair(cols[k],values[k]));
343 }
344
345 std::sort(work.begin(),work.end(),comp_pair_c);
346
347 for(unsigned int i=0 ; i<work.size() ; i++) {
348 ind.push_back(work[i].first);
349 x.push_back(work[i].second);
350 }
351 }
352
353 }
354
355 }
356
357
359 scmatrix(const srmatrix& A) : p(A.p), ind(A.ind), m(A.m), n(A.n), lb1(A.lb1), ub1(A.ub1), lb2(A.lb2), ub2(A.ub2) {
360 x.reserve(A.get_nnz());
361 for(unsigned int i=0 ; i<A.x.size() ; i++)
362 x.push_back(complex(A.x[i]));
363 }
364
365
367 scmatrix(const rmatrix& A) : m(ColLen(A)),n(RowLen(A)),lb1(Lb(A,1)),ub1(Ub(A,1)),lb2(Lb(A,2)),ub2(Ub(A,2)) {
368 p = std::vector<int>((n>0) ? n+1 : 1, 0);
369 ind.reserve((m*n*0.1 < 2*m) ? (int)(m*n*0.1) : 2*m);
370 x.reserve((m*n*0.1 < 2*m) ? (int)(m*n*0.1) : 2*m);
371
372 p[0] = 0;
373 int nnz = 0;
374
375 for(int j=0 ; j<n ; j++) {
376 for(int i=0 ; i<m ; i++) {
377 if(A[i+lb1][j+lb2] != 0.0) {
378 ind.push_back(i);
379 x.push_back(complex(A[i+lb1][j+lb2]));
380 nnz++;
381 }
382 }
383
384 p[j+1] = nnz;
385 }
386
387 }
388
390 scmatrix(const cmatrix& A) : m(ColLen(A)),n(RowLen(A)),lb1(Lb(A,1)),ub1(Ub(A,1)),lb2(Lb(A,2)),ub2(Ub(A,2)) {
391 p = std::vector<int>((n>0) ? n+1 : 1, 0);
392 ind.reserve((m*n*0.1 < 2*m) ? (int)(m*n*0.1) : 2*m);
393 x.reserve((m*n*0.1 < 2*m) ? (int)(m*n*0.1) : 2*m);
394
395 p[0] = 0;
396 int nnz = 0;
397
398 for(int j=0 ; j<n ; j++) {
399 for(int i=0 ; i<m ; i++) {
400 if(A[i+lb1][j+lb2] != 0.0) {
401 ind.push_back(i);
402 x.push_back(complex(A[i+lb1][j+lb2]));
403 nnz++;
404 }
405 }
406
407 p[j+1] = nnz;
408 }
409
410 }
411
413
416 scmatrix(const int ms, const int ns, const cmatrix& A) : m(ms), n(ns), lb1(1), ub1(ms), lb2(1), ub2(ns) {
417 //Banded matrix constructor
418 int nnz = RowLen(A)*ColLen(A);
419 p = std::vector<int>((n>0) ? n+1 : 1, 0);
420 ind.reserve(nnz);
421 x.reserve(nnz);
422
423 std::vector<triplet_store<complex> > work;
424 work.reserve(nnz);
425
426
427 for(int i=0 ; i<ColLen(A) ; i++) {
428 for(int j=Lb(A,2) ; j<=Ub(A,2) ; j++) {
429 if(i+j >=0 && i+j < n) {
430 work.push_back(triplet_store<complex>(i,i+j,A[i+Lb(A,1)][j]));
431 }
432 }
433 }
434
435 sort(work.begin(), work.end());
436
437 int i=0;
438
439 for(int j=0 ; j<n ; j++) {
440
441 while((unsigned int)i < work.size() && work[i].col == j ) {
442 ind.push_back(work[i].row);
443 x.push_back(work[i].val);
444 i++;
445 }
446
447 p[j+1] = i;
448 }
449
450 }
451
453 scmatrix(const srmatrix_slice&);
455 scmatrix(const scmatrix_slice&);
456
458 void full(cmatrix& A) const {
459 A = cmatrix(lb1,ub1,lb2,ub2);
460 A = 0.0;
461 for(int j=0 ; j<n ; j++) {
462 for(int k=p[j] ; k<p[j+1] ; k++) {
463 A[ind[k]+lb1][j+lb2] = x[k];
464 }
465 }
466 }
467
469
473 void dropzeros() {
474 std::vector<int> pnew(n+1,0);
475 std::vector<int> indnew;
476 std::vector<complex> xnew;
477 int nnznew = 0;
478
479 for(int j=0 ; j<n ; j++) {
480 for(int k=p[j] ; k<p[j+1] ; k++) {
481 if(x[k] != 0.0) {
482 xnew.push_back(x[k]);
483 indnew.push_back(ind[k]);
484 nnznew++;
485 }
486 }
487 pnew[j+1] = nnznew;
488 }
489
490 p = pnew;
491 ind = indnew;
492 x = xnew;
493 }
494
495
498 return sp_ms_assign<scmatrix,real,complex>(*this,A);
499 }
500
503 return sp_ms_assign<scmatrix,complex,complex>(*this,A);
504 }
505
508 return spf_mm_assign<scmatrix,rmatrix,complex>(*this,A);
509 }
510
513 return spf_mm_assign<scmatrix,cmatrix,complex>(*this,A);
514 }
515
518 return spf_mm_assign<scmatrix,rmatrix_slice,complex>(*this,A);
519 }
520
523 return spf_mm_assign<scmatrix,cmatrix_slice,complex>(*this,A);
524 }
525
528 m = A.m;
529 n = A.n;
530 p = A.p;
531 ind = A.ind;
532 x.clear();
533 x.reserve(A.get_nnz());
534 for(unsigned int i=0 ; i<A.x.size() ; i++)
535 x.push_back(complex(A.x[i]));
536 return *this;
537 }
538
539 /* scmatrix& operator=(const scmatrix& A) {
540 p = A.p;
541 ind = A.ind;
542 x = A.x;
543 return *this;
544 } */
545
550
552
558 const complex operator()(int i, int j) const {
559#if(CXSC_INDEX_CHECK)
560 if(i<lb1 || i>ub1 || j<lb2 || j>ub2)
561 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix::operator()(int, int)"));
562#endif
563 complex r(0.0);
564 for(int k=p[j-lb2] ; k<p[j-lb2+1] && ind[k]<=i-lb1 ; k++) {
565 if(ind[k] == i-lb1) r = x[k];
566 }
567 return r;
568 }
569
571
579 complex& element(int i, int j) {
580#if(CXSC_INDEX_CHECK)
581 if(i<lb1 || i>ub1 || j<lb2 || j>ub2)
582 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix::element()(int, int)"));
583#endif
584 int k;
585 for(k=p[j-lb2] ; k<p[j-lb2+1] && ind[k]<=i-lb1 ; k++) {
586 if(ind[k] == i-lb1) return x[k];
587 }
588
589 //Nicht gefunden, Element muss angelegt werden, da Schreibzugriff moeglich
590 std::vector<int>::iterator ind_it = ind.begin() + k;
591 std::vector<complex>::iterator x_it = x.begin() + k;
592 ind.insert(ind_it, i-lb1);
593 x_it = x.insert(x_it, complex(0.0));
594 for(k=j-lb2+1 ; k<(int)p.size() ; k++)
595 p[k]++;
596
597 return *x_it;
598 }
599
601 scmatrix_subv operator[](const cxscmatrix_column&);
603 scmatrix_subv operator[](const int);
605 const scmatrix_subv operator[](const cxscmatrix_column&) const;
607 const scmatrix_subv operator[](const int) const;
608
610 scmatrix_slice operator()(const int, const int , const int, const int);
612 const scmatrix_slice operator()(const int, const int , const int, const int) const;
613
615 scmatrix operator()(const intvector& pervec, const intvector& q) {
616 scmatrix A(m,n,get_nnz());
617 intvector per = perminv(pervec);
618
619 int nnz=0;
620 for(int k=0 ; k<n ; k++) {
621 A.p[k] = nnz;
622
623 std::map<int,complex> work;
624 for(int j=p[q[Lb(q)+k]] ; j<p[q[Lb(q)+k]+1] ; j++)
625 work.insert(std::make_pair(per[Lb(per)+ind[j]], x[j]));
626
627 for(std::map<int,complex>::iterator it = work.begin() ; it != work.end() ; it++) {
628 A.ind.push_back(it->first);
629 A.x.push_back(it->second);
630 }
631
632 nnz += work.size();
633
634 }
635
636 A.p[n] = nnz;
637
638 return A;
639 }
640
643 scmatrix A(m,n,get_nnz());
644 intvector per = perminv(pervec);
645
646 for(int k=0 ; k<n ; k++) {
647 A.p[k] = p[k];
648
649 std::map<int,complex> work;
650 for(int j=p[k] ; j<p[k+1] ; j++)
651 work.insert(std::make_pair(per[Lb(per)+ind[j]], x[j]));
652
653 for(std::map<int,complex>::iterator it = work.begin() ; it != work.end() ; it++) {
654 A.ind.push_back(it->first);
655 A.x.push_back(it->second);
656 }
657
658 }
659
660 A.p[n] = p[n];
661
662 return A;
663 }
664
667 intvector p = permvec(P);
668 intvector q = perminv(permvec(Q));
669 return (*this)(p,q);
670 }
671
674 intvector p = permvec(P);
675 return (*this)(p);
676 }
677
679 real density() const {
680 return p[n]/((double)m*n);
681 }
682
684 int get_nnz() const {
685 return p[n];
686 }
687
690 return spf_mm_addassign<scmatrix,rmatrix,cmatrix>(*this,B);
691 }
692
695 return spf_mm_addassign<scmatrix,cmatrix,cmatrix>(*this,B);
696 }
697
700 return spf_mm_addassign<scmatrix,rmatrix_slice,cmatrix>(*this,B);
701 }
702
705 return spf_mm_addassign<scmatrix,cmatrix_slice,cmatrix>(*this,B);
706 }
707
710 return spsp_mm_addassign<scmatrix,srmatrix,complex>(*this,B);
711 }
712
715 return spsp_mm_addassign<scmatrix,scmatrix,complex>(*this,B);
716 }
717
720 return spf_mm_subassign<scmatrix,rmatrix,cmatrix>(*this,B);
721 }
722
725 return spf_mm_subassign<scmatrix,cmatrix,cmatrix>(*this,B);
726 }
727
730 return spf_mm_subassign<scmatrix,rmatrix_slice,cmatrix>(*this,B);
731 }
732
735 return spf_mm_subassign<scmatrix,cmatrix_slice,cmatrix>(*this,B);
736 }
737
740 return spsp_mm_subassign<scmatrix,srmatrix,complex>(*this,B);
741 }
742
745 return spsp_mm_subassign<scmatrix,scmatrix,complex>(*this,B);
746 }
747
750 return spf_mm_multassign<scmatrix,cmatrix,sparse_cdot,cmatrix>(*this,B);
751 }
752
755 return spf_mm_multassign<scmatrix,rmatrix,sparse_cdot,cmatrix>(*this,B);
756 }
757
760 return spf_mm_multassign<scmatrix,rmatrix_slice,sparse_cdot,cmatrix>(*this,B);
761 }
762
765 return spf_mm_multassign<scmatrix,cmatrix_slice,sparse_cdot,cmatrix>(*this,B);
766 }
767
770 return spsp_mm_multassign<scmatrix,srmatrix,sparse_cdot,complex>(*this,B);
771 }
772
775 return spsp_mm_multassign<scmatrix,scmatrix,sparse_cdot,complex>(*this,B);
776 }
777
780 return sp_ms_multassign(*this,r);
781 }
782
785 return sp_ms_multassign(*this,r);
786 }
787
790 return sp_ms_divassign(*this,r);
791 }
792
795 return sp_ms_divassign(*this,r);
796 }
797
798 friend void SetLb(scmatrix&, const int, const int);
799 friend void SetUb(scmatrix&, const int, const int);
800 friend int Lb(const scmatrix&, int);
801 friend int Ub(const scmatrix&, int);
802 friend int RowLen(const scmatrix&);
803 friend int ColLen(const scmatrix&);
804 friend srmatrix Re(const scmatrix&);
805 friend srmatrix Im(const scmatrix&);
806 friend scmatrix Inf(const scimatrix&);
807 friend scmatrix Sup(const scimatrix&);
808 friend scmatrix mid(const scimatrix&);
809 friend scmatrix diam(const scimatrix&);
810 friend srmatrix abs(const scmatrix&);
811
812 friend srmatrix CompMat(const scmatrix&);
813 friend scmatrix transp(const scmatrix&);
814 friend scmatrix Id(const scmatrix&);
815
816 friend std::istream& operator>>(std::istream&, scmatrix_slice&);
817 friend std::istream& operator>>(std::istream&, scmatrix_subv&);
818
819 friend class srmatrix_slice;
820 friend class srmatrix_subv;
821 friend class srvector;
822 friend class scmatrix_slice;
823 friend class scmatrix_subv;
824 friend class scvector;
825 friend class scivector;
826 friend class scimatrix;
827 friend class scimatrix_slice;
828 friend class scimatrix_subv;
829 friend class cmatrix;
830 friend class cimatrix;
831
832
833#include "matrix_friend_declarations.inl"
834};
835
836inline cmatrix::cmatrix(const srmatrix& A) {
837 dat = new complex[A.m*A.n];
838 lb1 = A.lb1; lb2 = A.lb2; ub1 = A.ub1; ub2 = A.ub2;
839 xsize = A.n;
840 ysize = A.m;
841 *this = 0.0;
842 for(int j=0 ; j<A.n ; j++) {
843 for(int k=A.p[j] ; k<A.p[j+1] ; k++) {
844 dat[A.ind[k]*A.n+j] = A.x[k];
845 }
846 }
847}
848
849inline cmatrix::cmatrix(const scmatrix& A) {
850 dat = new complex[A.m*A.n];
851 lb1 = A.lb1; lb2 = A.lb2; ub1 = A.ub1; ub2 = A.ub2;
852 xsize = A.n;
853 ysize = A.m;
854 *this = 0.0;
855 for(int j=0 ; j<A.n ; j++) {
856 for(int k=A.p[j] ; k<A.p[j+1] ; k++) {
857 dat[A.ind[k]*A.n+j] = A.x[k];
858 }
859 }
860}
861
863inline scmatrix Id(const scmatrix& A) {
864 scmatrix I(A.m, A.n, (A.m>A.n) ? A.m : A.n);
865 I.lb1 = A.lb1; I.lb2 = A.lb2;
866 I.ub1 = A.ub1; I.ub2 = A.ub2;
867
868 if(A.m < A.n) {
869 for(int i=0 ; i<A.m ; i++) {
870 I.p[i+1] = I.p[i] + 1;
871 I.ind.push_back(i);
872 I.x.push_back(complex(1.0));
873 }
874 } else {
875 for(int i=0 ; i<A.n ; i++) {
876 I.p[i+1] = I.p[i] + 1;
877 I.ind.push_back(i);
878 I.x.push_back(complex(1.0));
879 }
880 }
881
882 return I;
883}
884
886inline scmatrix transp(const scmatrix& A) {
887 scmatrix B(A.n, A.m, A.get_nnz());
888
889 //Nichtnullen pro Zeile bestimmen
890 std::vector<int> w(A.m,0);
891 for(unsigned int i=0 ; i<A.ind.size() ; i++)
892 w[A.ind[i]]++;
893
894 //Spalten"pointer" setzen
895 B.p.resize(A.m+1);
896 B.p[0] = 0;
897 for(unsigned int i=1 ; i<B.p.size() ; i++)
898 B.p[i] = w[i-1] + B.p[i-1];
899
900 //w vorbereiten
901 w.insert(w.begin(), 0);
902 for(unsigned int i=1 ; i<w.size() ; i++) {
903 w[i] += w[i-1];
904 }
905
906 //neuer zeilenindex und wert wird gesetzt
907 int q;
908 B.ind.resize(A.get_nnz());
909 B.x.resize(A.get_nnz());
910 for(int j=0 ; j<A.n ; j++) {
911 for(int k=A.p[j] ; k<A.p[j+1] ; k++) {
912 q = w[A.ind[k]]++;
913 B.ind[q] = j;
914 B.x[q] = A.x[k];
915 }
916 }
917
918 return B;
919}
920
922
926inline void SetLb(scmatrix& A, const int i, const int j) {
927 if(i==1) {
928 A.lb1 = j;
929 A.ub1 = j + A.m - 1;
930 } else if(i==2) {
931 A.lb2 = j;
932 A.ub2 = j + A.n - 1;
933 }
934}
935
937
941inline void SetUb(scmatrix& A, const int i, const int j) {
942 if(i==1) {
943 A.ub1 = j;
944 A.lb1 = j - A.m + 1;
945 } else if(i==2) {
946 A.ub2 = j;
947 A.lb2 = j - A.n + 1;
948 }
949}
950
951
953
956inline int Lb(const scmatrix& A, int i) {
957 if(i==1)
958 return A.lb1;
959 else if(i==2)
960 return A.lb2;
961 else
962 return 1;
963}
964
966
969inline int Ub(const scmatrix& A, int i) {
970 if(i==1)
971 return A.ub1;
972 else if(i==2)
973 return A.ub2;
974 else
975 return 1;
976}
977
979inline int RowLen(const scmatrix& A) {
980 return A.n;
981}
982
984inline int ColLen(const scmatrix& A) {
985 return A.m;
986}
987
989inline void Resize(scmatrix& A) {
990 sp_m_resize(A);
991}
992
994inline void Resize(scmatrix& A, const int m, const int n) {
995 sp_m_resize(A,m,n);
996}
997
999inline void Resize(scmatrix& A, const int l1, const int u1, const int l2, const int u2) {
1000 sp_m_resize(A,l1,u1,l2,u2);
1001}
1002
1004inline srmatrix Re(const scmatrix& A) {
1005 srmatrix res(A.m,A.n,A.get_nnz());
1006 res.lb1 = A.lb1;
1007 res.lb2 = A.lb2;
1008 res.ub1 = A.ub1;
1009 res.ub2 = A.ub2;
1010 res.p = A.p;
1011 res.ind = A.ind;
1012
1013 for(int i=0 ; i<res.get_nnz() ; i++)
1014 res.x.push_back(Re(A.x[i]));
1015
1016 res.dropzeros();
1017
1018 return res;
1019}
1020
1022inline srmatrix Im(const scmatrix& A) {
1023 srmatrix res(A.m,A.n,A.get_nnz());
1024 res.lb1 = A.lb1;
1025 res.lb2 = A.lb2;
1026 res.ub1 = A.ub1;
1027 res.ub2 = A.ub2;
1028 res.p = A.p;
1029 res.ind = A.ind;
1030
1031 for(int i=0 ; i<res.get_nnz() ; i++)
1032 res.x.push_back(Im(A.x[i]));
1033
1034 res.dropzeros();
1035
1036 return res;
1037}
1038
1040inline srmatrix abs(const scmatrix& A) {
1041 srmatrix ret;
1042 ret.ind = A.ind;
1043 ret.p = A.p;
1044 for(unsigned int i=0 ; i<ret.x.size() ; i++)
1045 ret.x[i] = abs(A.x[i]);
1046 return ret;
1047}
1048
1050inline srmatrix CompMat(const scmatrix& A) {
1051 srmatrix res(A.m,A.n,A.get_nnz());
1052 res.lb1 = A.lb1;
1053 res.lb2 = A.lb2;
1054 res.ub1 = A.ub1;
1055 res.ub2 = A.ub2;
1056 res.p = A.p;
1057 res.ind = A.ind;
1058 res.p[A.n] = A.p[A.n];
1059
1060 for(int j=0 ; j<res.n ; j++) {
1061 for(int k=A.p[j] ; k<A.p[j+1] ; k++) {
1062 if(A.ind[k] == j)
1063 res.x.push_back(abs(A.x[k]));
1064 else
1065 res.x.push_back(-abs(A.x[k]));
1066 }
1067 }
1068
1069 res.dropzeros();
1070
1071 return res;
1072}
1073
1075
1081inline cmatrix operator*(const cmatrix& A, const srmatrix& B) {
1082 return fsp_mm_mult<cmatrix,srmatrix,cmatrix,sparse_cdot>(A,B);
1083}
1084
1086
1092inline cmatrix operator*(const rmatrix& A, const scmatrix& B) {
1093 return fsp_mm_mult<rmatrix,scmatrix,cmatrix,sparse_cdot>(A,B);
1094}
1095
1097
1103inline cmatrix operator*(const cmatrix& A, const scmatrix& B) {
1104 return fsp_mm_mult<cmatrix,scmatrix,cmatrix,sparse_cdot>(A,B);
1105}
1106
1108
1114inline cmatrix operator*(const scmatrix& A, const rmatrix& B) {
1115 return spf_mm_mult<scmatrix,rmatrix,cmatrix,sparse_cdot>(A,B);
1116}
1117
1119
1125inline cmatrix operator*(const srmatrix& A, const cmatrix& B) {
1126 return spf_mm_mult<srmatrix,cmatrix,cmatrix,sparse_cdot>(A,B);
1127}
1128
1130
1136inline cmatrix operator*(const scmatrix& A, const cmatrix& B) {
1137 return spf_mm_mult<scmatrix,cmatrix,cmatrix,sparse_cdot>(A,B);
1138}
1139
1141
1147inline cmatrix operator*(const cmatrix_slice& A, const srmatrix& B) {
1148 return fsp_mm_mult<cmatrix_slice,srmatrix,cmatrix,sparse_cdot>(A,B);
1149}
1150
1152
1158inline cmatrix operator*(const rmatrix_slice& A, const scmatrix& B) {
1159 return fsp_mm_mult<rmatrix_slice,scmatrix,cmatrix,sparse_cdot>(A,B);
1160}
1161
1163
1169inline cmatrix operator*(const cmatrix_slice& A, const scmatrix& B) {
1170 return fsp_mm_mult<cmatrix_slice,scmatrix,cmatrix,sparse_cdot>(A,B);
1171}
1172
1174
1180inline cmatrix operator*(const scmatrix& A, const rmatrix_slice& B) {
1181 return spf_mm_mult<scmatrix,rmatrix_slice,cmatrix,sparse_cdot>(A,B);
1182}
1183
1185
1191inline cmatrix operator*(const srmatrix& A, const cmatrix_slice& B) {
1192 return spf_mm_mult<srmatrix,cmatrix_slice,cmatrix,sparse_cdot>(A,B);
1193}
1194
1196
1202inline cmatrix operator*(const scmatrix& A, const cmatrix_slice& B) {
1203 return spf_mm_mult<scmatrix,cmatrix_slice,cmatrix,sparse_cdot>(A,B);
1204}
1205
1207
1213inline scmatrix operator*(const scmatrix& A, const srmatrix& B) {
1214 return spsp_mm_mult<scmatrix,srmatrix,scmatrix,sparse_cdot,complex>(A,B);
1215}
1216
1218
1224inline scmatrix operator*(const srmatrix& A, const scmatrix& B) {
1225 return spsp_mm_mult<srmatrix,scmatrix,scmatrix,sparse_cdot,complex>(A,B);
1226}
1227
1229
1235inline scmatrix operator*(const scmatrix& A, const scmatrix& B) {
1236 return spsp_mm_mult<scmatrix,scmatrix,scmatrix,sparse_cdot,complex>(A,B);
1237}
1238
1240inline scmatrix operator/(const scmatrix& A, const real& r) {
1241 return sp_ms_div<scmatrix,real,scmatrix>(A,r);
1242}
1243
1245inline scmatrix operator/(const scmatrix& A, const complex& r) {
1246 return sp_ms_div<scmatrix,complex,scmatrix>(A,r);
1247}
1248
1250inline scmatrix operator/(const srmatrix& A, const complex& r) {
1251 return sp_ms_div<srmatrix,complex,scmatrix>(A,r);
1252}
1253
1255inline scmatrix operator*(const scmatrix& A, const real& r) {
1256 return sp_ms_mult<scmatrix,real,scmatrix>(A,r);
1257}
1258
1260inline scmatrix operator*(const scmatrix& A, const complex& r) {
1261 return sp_ms_mult<scmatrix,complex,scmatrix>(A,r);
1262}
1263
1265inline scmatrix operator*(const srmatrix& A, const complex& r) {
1266 return sp_ms_mult<srmatrix,complex,scmatrix>(A,r);
1267}
1268
1270inline scmatrix operator*(const real& r, const scmatrix& A) {
1271 return sp_sm_mult<real,scmatrix,scmatrix>(r,A);
1272}
1273
1275inline scmatrix operator*(const complex& r, const scmatrix& A) {
1276 return sp_sm_mult<complex,scmatrix,scmatrix>(r,A);
1277}
1278
1280inline scmatrix operator*(const complex& r, const srmatrix& A) {
1281 return sp_sm_mult<complex,srmatrix,scmatrix>(r,A);
1282}
1283
1285
1291inline cvector operator*(const scmatrix& A, const rvector& v) {
1292 return spf_mv_mult<scmatrix,rvector,cvector,sparse_cdot>(A,v);
1293}
1294
1296
1302inline cvector operator*(const srmatrix& A, const cvector& v) {
1303 return spf_mv_mult<srmatrix,cvector,cvector,sparse_cdot>(A,v);
1304}
1305
1307
1313inline cvector operator*(const scmatrix& A, const cvector& v) {
1314 return spf_mv_mult<scmatrix,cvector,cvector,sparse_cdot>(A,v);
1315}
1316
1318
1324inline cvector operator*(const scmatrix& A, const rvector_slice& v) {
1325 return spf_mv_mult<scmatrix,rvector_slice,cvector,sparse_cdot>(A,v);
1326}
1327
1329
1335inline cvector operator*(const srmatrix& A, const cvector_slice& v) {
1336 return spf_mv_mult<srmatrix,cvector_slice,cvector,sparse_cdot>(A,v);
1337}
1338
1340
1346inline cvector operator*(const scmatrix& A, const cvector_slice& v) {
1347 return spf_mv_mult<scmatrix,cvector_slice,cvector,sparse_cdot>(A,v);
1348}
1349
1351
1357inline scvector operator*(const scmatrix& A, const srvector& v) {
1358 return spsp_mv_mult<scmatrix,srvector,scvector,sparse_cdot,complex>(A,v);
1359}
1360
1362
1368inline scvector operator*(const srmatrix& A, const scvector& v) {
1369 return spsp_mv_mult<srmatrix,scvector,scvector,sparse_cdot,complex>(A,v);
1370}
1371
1373
1379inline scvector operator*(const scmatrix& A, const scvector& v) {
1380 return spsp_mv_mult<scmatrix,scvector,scvector,sparse_cdot,complex>(A,v);
1381}
1382
1384
1390inline scvector operator*(const scmatrix& A, const srvector_slice& v) {
1391 return spsl_mv_mult<scmatrix,srvector_slice,scvector,sparse_cdot,complex>(A,v);
1392}
1393
1395
1401inline scvector operator*(const srmatrix& A, const scvector_slice& v) {
1402 return spsl_mv_mult<srmatrix,scvector_slice,scvector,sparse_cdot,complex>(A,v);
1403}
1404
1406
1412inline scvector operator*(const scmatrix& A, const scvector_slice& v) {
1413 return spsl_mv_mult<scmatrix,scvector_slice,scvector,sparse_cdot,complex>(A,v);
1414}
1415
1417
1423inline cvector operator*(const cmatrix& A, const srvector& v) {
1424 return fsp_mv_mult<cmatrix,srvector,cvector,sparse_cdot>(A,v);
1425}
1426
1428
1434inline cvector operator*(const rmatrix& A, const scvector& v) {
1435 return fsp_mv_mult<rmatrix,scvector,cvector,sparse_cdot>(A,v);
1436}
1437
1439
1445inline cvector operator*(const cmatrix& A, const scvector& v) {
1446 return fsp_mv_mult<cmatrix,scvector,cvector,sparse_cdot>(A,v);
1447}
1448
1450
1456inline cvector operator*(const cmatrix_slice& A, const srvector& v) {
1457 return fsp_mv_mult<cmatrix_slice,srvector,cvector,sparse_cdot>(A,v);
1458}
1459
1461
1467inline cvector operator*(const rmatrix_slice& A, const scvector& v) {
1468 return fsp_mv_mult<rmatrix_slice,scvector,cvector,sparse_cdot>(A,v);
1469}
1470
1472
1478inline cvector operator*(const cmatrix_slice& A, const scvector& v) {
1479 return fsp_mv_mult<cmatrix_slice,scvector,cvector,sparse_cdot>(A,v);
1480}
1481
1483
1489inline cvector operator*(const cmatrix& A, const srvector_slice& v) {
1490 return fsl_mv_mult<cmatrix,srvector_slice,cvector,sparse_cdot>(A,v);
1491}
1492
1494
1500inline cvector operator*(const rmatrix& A, const scvector_slice& v) {
1501 return fsl_mv_mult<rmatrix,scvector_slice,cvector,sparse_cdot>(A,v);
1502}
1503
1505
1511inline cvector operator*(const cmatrix& A, const scvector_slice& v) {
1512 return fsl_mv_mult<cmatrix,scvector_slice,cvector,sparse_cdot>(A,v);
1513}
1514
1516
1522inline cvector operator*(const cmatrix_slice& A, const srvector_slice& v) {
1523 return fsl_mv_mult<cmatrix_slice,srvector_slice,cvector,sparse_cdot>(A,v);
1524}
1525
1527
1533inline cvector operator*(const rmatrix_slice& A, const scvector_slice& v) {
1534 return fsl_mv_mult<rmatrix_slice,scvector_slice,cvector,sparse_cdot>(A,v);
1535}
1536
1538
1544inline cvector operator*(const cmatrix_slice& A, const scvector_slice& v) {
1545 return fsl_mv_mult<cmatrix_slice,scvector_slice,cvector,sparse_cdot>(A,v);
1546}
1547
1549inline cmatrix operator+(const cmatrix& A, const srmatrix& B) {
1550 return fsp_mm_add<cmatrix,srmatrix,cmatrix>(A,B);
1551}
1552
1554inline cmatrix operator+(const rmatrix& A, const scmatrix& B) {
1555 return fsp_mm_add<rmatrix,scmatrix,cmatrix>(A,B);
1556}
1557
1559inline cmatrix operator+(const cmatrix& A, const scmatrix& B) {
1560 return fsp_mm_add<cmatrix,scmatrix,cmatrix>(A,B);
1561}
1562
1564inline cmatrix operator+(const scmatrix& A, const rmatrix& B) {
1565 return spf_mm_add<scmatrix,rmatrix,cmatrix>(A,B);
1566}
1567
1569inline cmatrix operator+(const srmatrix& A, const cmatrix& B) {
1570 return spf_mm_add<srmatrix,cmatrix,cmatrix>(A,B);
1571}
1572
1574inline cmatrix operator+(const scmatrix& A, const cmatrix& B) {
1575 return spf_mm_add<scmatrix,cmatrix,cmatrix>(A,B);
1576}
1577
1579inline cmatrix operator+(const cmatrix_slice& A, const srmatrix& B) {
1580 return fsp_mm_add<cmatrix_slice,srmatrix,cmatrix>(A,B);
1581}
1582
1584inline cmatrix operator+(const rmatrix_slice& A, const scmatrix& B) {
1585 return fsp_mm_add<rmatrix_slice,scmatrix,cmatrix>(A,B);
1586}
1587
1589inline cmatrix operator+(const cmatrix_slice& A, const scmatrix& B) {
1590 return fsp_mm_add<cmatrix_slice,scmatrix,cmatrix>(A,B);
1591}
1592
1594inline cmatrix operator+(const scmatrix& A, const rmatrix_slice& B) {
1595 return spf_mm_add<scmatrix,rmatrix_slice,cmatrix>(A,B);
1596}
1597
1599inline cmatrix operator+(const srmatrix& A, const cmatrix_slice& B) {
1600 return spf_mm_add<srmatrix,cmatrix_slice,cmatrix>(A,B);
1601}
1602
1604inline cmatrix operator+(const scmatrix& A, const cmatrix_slice& B) {
1605 return spf_mm_add<scmatrix,cmatrix_slice,cmatrix>(A,B);
1606}
1607
1609inline scmatrix operator+(const scmatrix& A, const srmatrix& B) {
1610 return spsp_mm_add<scmatrix,srmatrix,scmatrix,complex>(A,B);
1611}
1612
1614inline scmatrix operator+(const srmatrix& A, const scmatrix& B) {
1615 return spsp_mm_add<srmatrix,scmatrix,scmatrix,complex>(A,B);
1616}
1617
1619inline scmatrix operator+(const scmatrix& A, const scmatrix& B) {
1620 return spsp_mm_add<scmatrix,scmatrix,scmatrix,complex>(A,B);
1621}
1622
1624inline cmatrix operator-(const cmatrix& A, const srmatrix& B) {
1625 return fsp_mm_sub<cmatrix,srmatrix,cmatrix>(A,B);
1626}
1627
1629inline cmatrix operator-(const rmatrix& A, const scmatrix& B) {
1630 return fsp_mm_sub<rmatrix,scmatrix,cmatrix>(A,B);
1631}
1632
1634inline cmatrix operator-(const cmatrix& A, const scmatrix& B) {
1635 return fsp_mm_sub<cmatrix,scmatrix,cmatrix>(A,B);
1636}
1637
1639inline cmatrix operator-(const scmatrix& A, const rmatrix& B) {
1640 return spf_mm_sub<scmatrix,rmatrix,cmatrix>(A,B);
1641}
1642
1644inline cmatrix operator-(const srmatrix& A, const cmatrix& B) {
1645 return spf_mm_sub<srmatrix,cmatrix,cmatrix>(A,B);
1646}
1647
1649inline cmatrix operator-(const scmatrix& A, const cmatrix& B) {
1650 return spf_mm_sub<scmatrix,cmatrix,cmatrix>(A,B);
1651}
1652
1654inline cmatrix operator-(const cmatrix_slice& A, const srmatrix& B) {
1655 return fsp_mm_sub<cmatrix_slice,srmatrix,cmatrix>(A,B);
1656}
1657
1659inline cmatrix operator-(const rmatrix_slice& A, const scmatrix& B) {
1660 return fsp_mm_sub<rmatrix_slice,scmatrix,cmatrix>(A,B);
1661}
1662
1664inline cmatrix operator-(const cmatrix_slice& A, const scmatrix& B) {
1665 return fsp_mm_sub<cmatrix_slice,scmatrix,cmatrix>(A,B);
1666}
1667
1669inline cmatrix operator-(const scmatrix& A, const rmatrix_slice& B) {
1670 return spf_mm_sub<scmatrix,rmatrix_slice,cmatrix>(A,B);
1671}
1672
1674inline cmatrix operator-(const srmatrix& A, const cmatrix_slice& B) {
1675 return spf_mm_sub<srmatrix,cmatrix_slice,cmatrix>(A,B);
1676}
1677
1679inline cmatrix operator-(const scmatrix& A, const cmatrix_slice& B) {
1680 return spf_mm_sub<scmatrix,cmatrix_slice,cmatrix>(A,B);
1681}
1682
1684inline scmatrix operator-(const scmatrix& A, const srmatrix& B) {
1685 return spsp_mm_sub<scmatrix,srmatrix,scmatrix,complex>(A,B);
1686}
1687
1689inline scmatrix operator-(const srmatrix& A, const scmatrix& B) {
1690 return spsp_mm_sub<srmatrix,scmatrix,scmatrix,complex>(A,B);
1691}
1692
1694inline scmatrix operator-(const scmatrix& A, const scmatrix& B) {
1695 return spsp_mm_sub<scmatrix,scmatrix,scmatrix,complex>(A,B);
1696}
1697
1699inline scmatrix operator-(const scmatrix& M) {
1700 return sp_m_negative<scmatrix,scmatrix>(M);
1701}
1702
1704inline scmatrix& operator+(scmatrix& A) {
1705 return A;
1706}
1707
1709 *this = rmatrix(B);
1710 return *this;
1711}
1712
1714 *this = rmatrix(B);
1715 return *this;
1716}
1717
1719 *this = cmatrix(B);
1720 return *this;
1721}
1722
1724 *this = cmatrix(B);
1725 return *this;
1726}
1727
1729 return fsp_mm_addassign(*this,B);
1730}
1731
1733 return fsp_mm_addassign(*this,B);
1734}
1735
1737 return fsp_mm_addassign(*this,B);
1738}
1739
1741 return fsp_mm_addassign(*this,B);
1742}
1743
1745 return fsp_mm_subassign(*this,B);
1746}
1747
1749 return fsp_mm_subassign(*this,B);
1750}
1751
1753 return fsp_mm_subassign(*this,B);
1754}
1755
1757 return fsp_mm_subassign(*this,B);
1758}
1759
1761 return fsp_mm_multassign<cmatrix,srmatrix,sparse_cdot,cmatrix>(*this,B);
1762}
1763
1765 return fsp_mm_multassign<cmatrix,scmatrix,sparse_cdot,cmatrix>(*this,B);
1766}
1767
1769 return fsp_mm_multassign<cmatrix_slice,srmatrix,sparse_cdot,cmatrix>(*this,B);
1770}
1771
1773 return fsp_mm_multassign<cmatrix_slice,scmatrix,sparse_cdot,cmatrix>(*this,B);
1774}
1775
1777inline bool operator==(const scmatrix& A, const srmatrix& B) {
1778 return spsp_mm_comp(A,B);
1779}
1780
1782inline bool operator==(const srmatrix& A, const scmatrix& B) {
1783 return spsp_mm_comp(A,B);
1784}
1785
1787inline bool operator==(const scmatrix& A, const scmatrix& B) {
1788 return spsp_mm_comp(A,B);
1789}
1790
1792inline bool operator==(const scmatrix& A, const rmatrix& B) {
1793 return spf_mm_comp(A,B);
1794}
1795
1797inline bool operator==(const srmatrix& A, const cmatrix& B) {
1798 return spf_mm_comp(A,B);
1799}
1800
1802inline bool operator==(const scmatrix& A, const cmatrix& B) {
1803 return spf_mm_comp(A,B);
1804}
1805
1807inline bool operator==(const cmatrix& A, const srmatrix& B) {
1808 return fsp_mm_comp(A,B);
1809}
1810
1812inline bool operator==(const rmatrix& A, const scmatrix& B) {
1813 return fsp_mm_comp(A,B);
1814}
1815
1817inline bool operator==(const cmatrix& A, const scmatrix& B) {
1818 return fsp_mm_comp(A,B);
1819}
1820
1822inline bool operator==(const cmatrix_slice& A, const srmatrix& B) {
1823 return fsp_mm_comp(A,B);
1824}
1825
1827inline bool operator==(const rmatrix_slice& A, const scmatrix& B) {
1828 return fsp_mm_comp(A,B);
1829}
1830
1832inline bool operator==(const cmatrix_slice& A, const scmatrix& B) {
1833 return fsp_mm_comp(A,B);
1834}
1835
1837inline bool operator==(const scmatrix& A, const rmatrix_slice& B) {
1838 return spf_mm_comp(A,B);
1839}
1840
1842inline bool operator==(const srmatrix& A, const cmatrix_slice& B) {
1843 return spf_mm_comp(A,B);
1844}
1845
1847inline bool operator==(const scmatrix& A, const cmatrix_slice& B) {
1848 return spf_mm_comp(A,B);
1849}
1850
1852inline bool operator!=(const scmatrix& A, const srmatrix& B) {
1853 return !spsp_mm_comp(A,B);
1854}
1855
1857inline bool operator!=(const srmatrix& A, const scmatrix& B) {
1858 return !spsp_mm_comp(A,B);
1859}
1860
1862inline bool operator!=(const scmatrix& A, const scmatrix& B) {
1863 return !spsp_mm_comp(A,B);
1864}
1865
1867inline bool operator!=(const scmatrix& A, const rmatrix& B) {
1868 return !spf_mm_comp(A,B);
1869}
1870
1872inline bool operator!=(const srmatrix& A, const cmatrix& B) {
1873 return !spf_mm_comp(A,B);
1874}
1875
1877inline bool operator!=(const scmatrix& A, const cmatrix& B) {
1878 return !spf_mm_comp(A,B);
1879}
1880
1882inline bool operator!=(const cmatrix& A, const srmatrix& B) {
1883 return !fsp_mm_comp(A,B);
1884}
1885
1887inline bool operator!=(const rmatrix& A, const scmatrix& B) {
1888 return !fsp_mm_comp(A,B);
1889}
1890
1892inline bool operator!=(const cmatrix& A, const scmatrix& B) {
1893 return !fsp_mm_comp(A,B);
1894}
1895
1897inline bool operator!=(const cmatrix_slice& A, const srmatrix& B) {
1898 return !fsp_mm_comp(A,B);
1899}
1900
1902inline bool operator!=(const rmatrix_slice& A, const scmatrix& B) {
1903 return !fsp_mm_comp(A,B);
1904}
1905
1907inline bool operator!=(const cmatrix_slice& A, const scmatrix& B) {
1908 return !fsp_mm_comp(A,B);
1909}
1910
1912inline bool operator!=(const scmatrix& A, const rmatrix_slice& B) {
1913 return !spf_mm_comp(A,B);
1914}
1915
1917inline bool operator!=(const srmatrix& A, const cmatrix_slice& B) {
1918 return !spf_mm_comp(A,B);
1919}
1920
1922inline bool operator!=(const scmatrix& A, const cmatrix_slice& B) {
1923 return !spf_mm_comp(A,B);
1924}
1925
1927inline bool operator!(const scmatrix& A) {
1928 return sp_m_not(A);
1929}
1930
1932
1937inline std::ostream& operator<<(std::ostream& os, const scmatrix& A) {
1938 return sp_m_output<scmatrix,complex>(os,A);
1939}
1940
1942
1947inline std::istream& operator>>(std::istream& is, scmatrix& A) {
1948 return sp_m_input<scmatrix,complex>(is,A);
1949}
1950
1952
1957 public:
1958 scmatrix A;
1959 scmatrix* M; //Originalmatrix
1960
1961 private:
1962 scmatrix_slice(scmatrix& Mat, int sl1l, int sl1u, int sl2l, int sl2u) {
1963 A.lb1 = sl1l;
1964 A.lb2 = sl2l;
1965 A.ub1 = sl1u;
1966 A.ub2 = sl2u;
1967 A.m = sl1u-sl1l+1;
1968 A.n = sl2u-sl2l+1;
1969
1970 //Kopieren der Werte aus A
1971 A.p = std::vector<int>(A.n+1, 0);
1972 A.ind.reserve(A.m + A.n);
1973 A.x.reserve(A.m + A.n);
1974
1975 for(int i=0 ; i<A.n ; i++) {
1976 A.p[i+1] = A.p[i];
1977 for(int j=Mat.p[sl2l-Mat.lb2+i] ; j<Mat.p[sl2l-Mat.lb2+i+1] ; j++) {
1978 if(Mat.ind[j] >= sl1l-Mat.lb1 && Mat.ind[j] <= sl1u-Mat.lb1) {
1979 A.ind.push_back(Mat.ind[j]-(sl1l-Mat.lb1));
1980 A.x.push_back(Mat.x[j]);
1981 A.p[i+1]++;
1982 }
1983 }
1984 }
1985
1986 //Zeiger auf A fuer Datenmanipulationen
1987 M = &Mat;
1988 }
1989
1990 scmatrix_slice(const scmatrix& Mat, int sl1l, int sl1u, int sl2l, int sl2u) {
1991 A.lb1 = sl1l;
1992 A.lb2 = sl2l;
1993 A.ub1 = sl1u;
1994 A.ub2 = sl2u;
1995 A.m = sl1u-sl1l+1;
1996 A.n = sl2u-sl2l+1;
1997
1998 //Kopieren der Werte aus A
1999 A.p = std::vector<int>(A.n+1, 0);
2000 A.ind.reserve(A.m + A.n);
2001 A.x.reserve(A.m + A.n);
2002
2003 for(int i=0 ; i<A.n ; i++) {
2004 A.p[i+1] = A.p[i];
2005 for(int j=Mat.p[sl2l-Mat.lb2+i] ; j<Mat.p[sl2l-Mat.lb2+i+1] ; j++) {
2006 if(Mat.ind[j] >= sl1l-Mat.lb1 && Mat.ind[j] <= sl1u-Mat.lb1) {
2007 A.ind.push_back(Mat.ind[j]-(sl1l-Mat.lb1));
2008 A.x.push_back(Mat.x[j]);
2009 A.p[i+1]++;
2010 }
2011 }
2012 }
2013
2014 //Zeiger auf A fuer Datenmanipulationen
2015 M = const_cast<scmatrix*>(&Mat);
2016 }
2017
2018 public:
2021 return sl_ms_assign<scmatrix_slice, real, std::vector<complex>::iterator, complex>(*this,C);
2022 }
2023
2026 return sl_ms_assign<scmatrix_slice, complex, std::vector<complex>::iterator, complex>(*this,C);
2027 }
2028
2031 return slsp_mm_assign<scmatrix_slice, srmatrix, std::vector<complex>::iterator>(*this,C);
2032 }
2033
2036 return slsp_mm_assign<scmatrix_slice, scmatrix, std::vector<complex>::iterator>(*this,C);
2037 }
2038
2041 return slf_mm_assign<scmatrix_slice, rmatrix, std::vector<complex>::iterator, complex>(*this,C);
2042 }
2043
2046 return slf_mm_assign<scmatrix_slice, cmatrix, std::vector<complex>::iterator, complex>(*this,C);
2047 }
2048
2051 return slf_mm_assign<scmatrix_slice, rmatrix_slice, std::vector<complex>::iterator, complex>(*this,C);
2052 }
2053
2056 return slf_mm_assign<scmatrix_slice, cmatrix_slice, std::vector<complex>::iterator, complex>(*this,C);
2057 }
2058
2061 *this = C.A;
2062 return *this;
2063 }
2064
2067 *this = C.A;
2068 return *this;
2069 }
2070
2073 *this = A*M.A;
2074 return *this;
2075 }
2076
2079 *this = A*M.A;
2080 return *this;
2081 }
2082
2085 *this = A*M;
2086 return *this;
2087 }
2088
2091 *this = A*M;
2092 return *this;
2093 }
2094
2097 *this = A*M;
2098 return *this;
2099 }
2100
2103 *this = A*M;
2104 return *this;
2105 }
2106
2109 *this = A*M;
2110 return *this;
2111 }
2112
2115 *this = A*M;
2116 return *this;
2117 }
2118
2121 *this = A*r;
2122 return *this;
2123 }
2124
2127 *this = A*r;
2128 return *this;
2129 }
2130
2133 *this = A/r;
2134 return *this;
2135 }
2136
2139 *this = A/r;
2140 return *this;
2141 }
2142
2145 *this = A+M.A;
2146 return *this;
2147 }
2148
2151 *this = A+M.A;
2152 return *this;
2153 }
2154
2157 *this = A+M;
2158 return *this;
2159 }
2160
2163 *this = A+M;
2164 return *this;
2165 }
2166
2169 *this = A+M;
2170 return *this;
2171 }
2172
2175 *this = A+M;
2176 return *this;
2177 }
2178
2181 *this = A+M;
2182 return *this;
2183 }
2184
2187 *this = A+M;
2188 return *this;
2189 }
2190
2193 *this = A-M.A;
2194 return *this;
2195 }
2196
2199 *this = A-M.A;
2200 return *this;
2201 }
2202
2205 *this = A-M;
2206 return *this;
2207 }
2208
2211 *this = A-M;
2212 return *this;
2213 }
2214
2217 *this = A-M;
2218 return *this;
2219 }
2220
2223 *this = A-M;
2224 return *this;
2225 }
2226
2229 *this = A-M;
2230 return *this;
2231 }
2232
2235 *this = A-M;
2236 return *this;
2237 }
2238
2240
2244 const complex operator()(const int i, const int j) const {
2245#if(CXSC_INDEX_CHECK)
2246 if(i<A.lb1 || i>A.ub1 || j<A.lb2 || j>A.ub2)
2247 cxscthrow(ELEMENT_NOT_IN_VEC("scmatrix_slice::operator()(int, int)"));
2248#endif
2249 complex r = A(i,j);
2250 return r;
2251 }
2252
2254
2258 complex& element(const int i, const int j) {
2259#if(CXSC_INDEX_CHECK)
2260 if(i<A.lb1 || i>A.ub1 || j<A.lb2 || j>A.ub2)
2261 cxscthrow(ELEMENT_NOT_IN_VEC("scmatrix_slice::element(int, int)"));
2262#endif
2263 return M->element(i,j);
2264 }
2265
2267 scmatrix_subv operator[](const int);
2269 scmatrix_subv operator[](const cxscmatrix_column&);
2271 const scmatrix_subv operator[](const int) const;
2273 const scmatrix_subv operator[](const cxscmatrix_column&) const;
2274
2275 friend std::ostream& operator<<(std::ostream&, const scmatrix_slice&);
2276 friend std::istream& operator>>(std::istream&, const scmatrix_subv&);
2277
2278 friend int Lb(const scmatrix_slice&, const int);
2279 friend int Ub(const scmatrix_slice&, const int);
2280 friend srmatrix Re(const scmatrix_slice&);
2281 friend srmatrix Im(const scmatrix_slice&);
2282 friend int RowLen(const scmatrix_slice&);
2283 friend int ColLen(const scmatrix_slice&);
2284
2285 friend class srmatrix;
2286 friend class srmatrix_subv;
2287 friend class srvector;
2288 friend class scmatrix;
2289 friend class scmatrix_subv;
2290 friend class scvector;
2291 friend class scimatrix;
2292 friend class scimatrix_subv;
2293 friend class scimatrix_slice;
2294 friend class scivector;
2295 friend class cmatrix;
2296 friend class cimatrix;
2297
2298#include "matrix_friend_declarations.inl"
2299};
2300
2302 dat = new complex[A.A.m*A.A.n];
2303 lb1 = A.A.lb1; lb2 = A.A.lb2; ub1 = A.A.ub1; ub2 = A.A.ub2;
2304 xsize = A.A.n;
2305 ysize = A.A.m;
2306 *this = 0.0;
2307 for(int j=0 ; j<A.A.n ; j++) {
2308 for(int k=A.A.p[j] ; k<A.A.p[j+1] ; k++) {
2309 dat[A.A.ind[k]*A.A.n+j] = A.A.x[k];
2310 }
2311 }
2312}
2313
2315 dat = new complex[A.A.m*A.A.n];
2316 lb1 = A.A.lb1; lb2 = A.A.lb2; ub1 = A.A.ub1; ub2 = A.A.ub2;
2317 xsize = A.A.n;
2318 ysize = A.A.m;
2319 *this = 0.0;
2320 for(int j=0 ; j<A.A.n ; j++) {
2321 for(int k=A.A.p[j] ; k<A.A.p[j+1] ; k++) {
2322 dat[A.A.ind[k]*A.A.n+j] = A.A.x[k];
2323 }
2324 }
2325}
2326
2328inline int RowLen(const scmatrix_slice& S) {
2329 return RowLen(S.A);
2330}
2331
2333inline int ColLen(const scmatrix_slice& S) {
2334 return ColLen(S.A);
2335}
2336
2337inline scmatrix_slice scmatrix::operator()(const int i, const int j, const int k, const int l) {
2338#if(CXSC_INDEX_CHECK)
2339 if(i<lb1 || j>ub1 || k<lb2 || l>ub2)
2340 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix::operator()(int, int, int, int)"));
2341#endif
2342 return scmatrix_slice(*this, i, j, k, l);
2343}
2344
2345inline const scmatrix_slice scmatrix::operator()(const int i, const int j, const int k, const int l) const {
2346#if(CXSC_INDEX_CHECK)
2347 if(i<lb1 || j>ub1 || k<lb2 || l>ub2)
2348 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix::operator()(int, int, int, int) const"));
2349#endif
2350 return scmatrix_slice(*this, i, j, k, l);
2351}
2352
2354 *this = S.A;
2355 return *this;
2356}
2357
2359 *this = S.A;
2360 return *this;
2361}
2362
2364inline int Lb(const scmatrix_slice& S, const int i) {
2365 return Lb(S.A, i);
2366}
2367
2369inline int Ub(const scmatrix_slice& S, const int i) {
2370 return Ub(S.A, i);
2371}
2372
2374inline srmatrix Re(const scmatrix_slice& S) {
2375 return Re(S.A);
2376}
2377
2379inline srmatrix Im(const scmatrix_slice& S) {
2380 return Im(S.A);
2381}
2382
2384 m = S.A.m;
2385 n = S.A.n;
2386 lb1 = S.A.lb1;
2387 ub1 = S.A.ub1;
2388 lb2 = S.A.lb2;
2389 ub2 = S.A.ub2;
2390 *this = S.A;
2391}
2392
2394 m = S.A.m;
2395 n = S.A.n;
2396 lb1 = S.A.lb1;
2397 ub1 = S.A.ub1;
2398 lb2 = S.A.lb2;
2399 ub2 = S.A.ub2;
2400 *this = S.A;
2401}
2402
2404inline scmatrix operator-(const scmatrix_slice& M) {
2405 return sp_m_negative<scmatrix,scmatrix>(M.A);
2406}
2407
2409inline scmatrix operator+(const scmatrix_slice& M) {
2410 return M.A;
2411}
2412
2414
2420inline scmatrix operator*(const scmatrix_slice& M1, const srmatrix_slice& M2) {
2421 return spsp_mm_mult<scmatrix,srmatrix,scmatrix,sparse_cdot,complex>(M1.A,M2.A);
2422}
2423
2425
2431inline scmatrix operator*(const srmatrix_slice& M1, const scmatrix_slice& M2) {
2432 return spsp_mm_mult<srmatrix,scmatrix,scmatrix,sparse_cdot,complex>(M1.A,M2.A);
2433}
2434
2436
2442inline scmatrix operator*(const scmatrix_slice& M1, const scmatrix_slice& M2) {
2443 return spsp_mm_mult<scmatrix,scmatrix,scmatrix,sparse_cdot,complex>(M1.A,M2.A);
2444}
2445
2447
2453inline scmatrix operator*(const scmatrix_slice& M1, const srmatrix& M2) {
2454 return spsp_mm_mult<scmatrix,srmatrix,scmatrix,sparse_cdot,complex>(M1.A,M2);
2455}
2456
2458
2464inline scmatrix operator*(const srmatrix_slice& M1, const scmatrix& M2) {
2465 return spsp_mm_mult<srmatrix,scmatrix,scmatrix,sparse_cdot,complex>(M1.A,M2);
2466}
2467
2469
2475inline scmatrix operator*(const scmatrix_slice& M1, const scmatrix& M2) {
2476 return spsp_mm_mult<scmatrix,scmatrix,scmatrix,sparse_cdot,complex>(M1.A,M2);
2477}
2478
2480
2486inline scmatrix operator*(const scmatrix& M1, const srmatrix_slice& M2) {
2487 return spsp_mm_mult<scmatrix,srmatrix,scmatrix,sparse_cdot,complex>(M1,M2.A);
2488}
2489
2491
2497inline scmatrix operator*(const srmatrix& M1, const scmatrix_slice& M2) {
2498 return spsp_mm_mult<srmatrix,scmatrix,scmatrix,sparse_cdot,complex>(M1,M2.A);
2499}
2500
2502
2508inline scmatrix operator*(const scmatrix& M1, const scmatrix_slice& M2) {
2509 return spsp_mm_mult<scmatrix,scmatrix,scmatrix,sparse_cdot,complex>(M1,M2.A);
2510}
2511
2513
2519inline cmatrix operator*(const scmatrix_slice& M1, const rmatrix& M2) {
2520 return spf_mm_mult<scmatrix,rmatrix,cmatrix,sparse_cdot>(M1.A,M2);
2521}
2522
2524
2530inline cmatrix operator*(const srmatrix_slice& M1, const cmatrix& M2) {
2531 return spf_mm_mult<srmatrix,cmatrix,cmatrix,sparse_cdot>(M1.A,M2);
2532}
2533
2535
2541inline cmatrix operator*(const scmatrix_slice& M1, const cmatrix& M2) {
2542 return spf_mm_mult<scmatrix,cmatrix,cmatrix,sparse_cdot>(M1.A,M2);
2543}
2544
2546
2552inline cmatrix operator*(const cmatrix& M1, const srmatrix_slice& M2) {
2553 return fsp_mm_mult<cmatrix,srmatrix,cmatrix,sparse_cdot>(M1,M2.A);
2554}
2555
2557
2563inline cmatrix operator*(const rmatrix& M1, const scmatrix_slice& M2) {
2564 return fsp_mm_mult<rmatrix,scmatrix,cmatrix,sparse_cdot>(M1,M2.A);
2565}
2566
2568
2574inline cmatrix operator*(const cmatrix& M1, const scmatrix_slice& M2) {
2575 return fsp_mm_mult<cmatrix,scmatrix,cmatrix,sparse_cdot>(M1,M2.A);
2576}
2577
2579
2585inline cmatrix operator*(const scmatrix_slice& M1, const rmatrix_slice& M2) {
2586 return spf_mm_mult<scmatrix,rmatrix_slice,cmatrix,sparse_cdot>(M1.A,M2);
2587}
2588
2590
2596inline cmatrix operator*(const srmatrix_slice& M1, const cmatrix_slice& M2) {
2597 return spf_mm_mult<srmatrix,cmatrix_slice,cmatrix,sparse_cdot>(M1.A,M2);
2598}
2599
2601
2607inline cmatrix operator*(const scmatrix_slice& M1, const cmatrix_slice& M2) {
2608 return spf_mm_mult<scmatrix,cmatrix_slice,cmatrix,sparse_cdot>(M1.A,M2);
2609}
2610
2612
2618inline cmatrix operator*(const cmatrix_slice& M1, const srmatrix_slice& M2) {
2619 return fsp_mm_mult<cmatrix,srmatrix,cmatrix,sparse_cdot>(M1,M2.A);
2620}
2621
2623
2629inline cmatrix operator*(const rmatrix_slice& M1, const scmatrix_slice& M2) {
2630 return fsp_mm_mult<rmatrix,scmatrix,cmatrix,sparse_cdot>(M1,M2.A);
2631}
2632
2634
2640inline cmatrix operator*(const cmatrix_slice& M1, const scmatrix_slice& M2) {
2641 return fsp_mm_mult<cmatrix,scmatrix,cmatrix,sparse_cdot>(M1,M2.A);
2642}
2643
2645
2651inline scvector operator*(const scmatrix_slice& M, const srvector& v) {
2652 return spsp_mv_mult<scmatrix,srvector,scvector,sparse_cdot,complex>(M.A,v);
2653}
2654
2656
2662inline scvector operator*(const srmatrix_slice& M, const scvector& v) {
2663 return spsp_mv_mult<srmatrix,scvector,scvector,sparse_cdot,complex>(M.A,v);
2664}
2665
2667
2673inline scvector operator*(const scmatrix_slice& M, const scvector& v) {
2674 return spsp_mv_mult<scmatrix,scvector,scvector,sparse_cdot,complex>(M.A,v);
2675}
2676
2678
2685 return spsl_mv_mult<scmatrix,srvector_slice,scvector,sparse_cdot,complex>(M.A,v);
2686}
2687
2689
2696 return spsl_mv_mult<srmatrix,scvector_slice,scvector,sparse_cdot,complex>(M.A,v);
2697}
2698
2700
2707 return spsl_mv_mult<scmatrix,scvector_slice,scvector,sparse_cdot,complex>(M.A,v);
2708}
2709
2711
2717inline cvector operator*(const scmatrix_slice& M, const rvector& v) {
2718 return spf_mv_mult<scmatrix,rvector,cvector,sparse_cdot>(M.A,v);
2719}
2720
2722
2728inline cvector operator*(const srmatrix_slice& M, const cvector& v) {
2729 return spf_mv_mult<srmatrix,cvector,cvector,sparse_cdot>(M.A,v);
2730}
2731
2733
2739inline cvector operator*(const scmatrix_slice& M, const cvector& v) {
2740 return spf_mv_mult<scmatrix,cvector,cvector,sparse_cdot>(M.A,v);
2741}
2742
2744
2750inline cvector operator*(const scmatrix_slice& M, const rvector_slice& v) {
2751 return spf_mv_mult<scmatrix,rvector_slice,cvector,sparse_cdot>(M.A,v);
2752}
2753
2755
2761inline cvector operator*(const srmatrix_slice& M, const cvector_slice& v) {
2762 return spf_mv_mult<srmatrix,cvector_slice,cvector,sparse_cdot>(M.A,v);
2763}
2764
2766
2772inline cvector operator*(const scmatrix_slice& M, const cvector_slice& v) {
2773 return spf_mv_mult<scmatrix,cvector_slice,cvector,sparse_cdot>(M.A,v);
2774}
2775
2777inline scmatrix operator*(const scmatrix_slice& M, const real& r) {
2778 return sp_ms_mult<scmatrix,real,scmatrix>(M.A,r);
2779}
2780
2782inline scmatrix operator*(const scmatrix_slice& M, const complex& r) {
2783 return sp_ms_mult<scmatrix,complex,scmatrix>(M.A,r);
2784}
2785
2787inline scmatrix operator*(const srmatrix_slice& M, const complex& r) {
2788 return sp_ms_mult<srmatrix,complex,scmatrix>(M.A,r);
2789}
2790
2792inline scmatrix operator/(const scmatrix_slice& M, const real& r) {
2793 return sp_ms_div<scmatrix,real,scmatrix>(M.A,r);
2794}
2795
2797inline scmatrix operator/(const scmatrix_slice& M, const complex& r) {
2798 return sp_ms_div<scmatrix,complex,scmatrix>(M.A,r);
2799}
2800
2802inline scmatrix operator/(const srmatrix_slice& M, const complex& r) {
2803 return sp_ms_div<srmatrix,complex,scmatrix>(M.A,r);
2804}
2805
2807inline scmatrix operator*(const real& r, const scmatrix_slice& M) {
2808 return sp_sm_mult<real,scmatrix,scmatrix>(r,M.A);
2809}
2810
2812inline scmatrix operator*(const complex& r, const srmatrix_slice& M) {
2813 return sp_sm_mult<complex,srmatrix,scmatrix>(r,M.A);
2814}
2815
2817inline scmatrix operator*(const complex& r, const scmatrix_slice& M) {
2818 return sp_sm_mult<complex,scmatrix,scmatrix>(r,M.A);
2819}
2820
2822inline scmatrix operator+(const scmatrix_slice& M1, const srmatrix_slice& M2) {
2823 return spsp_mm_add<scmatrix,srmatrix,scmatrix,complex>(M1.A,M2.A);
2824}
2825
2827inline scmatrix operator+(const srmatrix_slice& M1, const scmatrix_slice& M2) {
2828 return spsp_mm_add<srmatrix,scmatrix,scmatrix,complex>(M1.A,M2.A);
2829}
2830
2832inline scmatrix operator+(const scmatrix_slice& M1, const scmatrix_slice& M2) {
2833 return spsp_mm_add<scmatrix,scmatrix,scmatrix,complex>(M1.A,M2.A);
2834}
2835
2837inline scmatrix operator+(const scmatrix_slice& M1, const srmatrix& M2) {
2838 return spsp_mm_add<scmatrix,srmatrix,scmatrix,complex>(M1.A,M2);
2839}
2840
2842inline scmatrix operator+(const srmatrix_slice& M1, const scmatrix& M2) {
2843 return spsp_mm_add<srmatrix,scmatrix,scmatrix,complex>(M1.A,M2);
2844}
2845
2847inline scmatrix operator+(const scmatrix_slice& M1, const scmatrix& M2) {
2848 return spsp_mm_add<scmatrix,scmatrix,scmatrix,complex>(M1.A,M2);
2849}
2850
2852inline scmatrix operator+(const scmatrix& M1, const srmatrix_slice& M2) {
2853 return spsp_mm_add<scmatrix,srmatrix,scmatrix,complex>(M1,M2.A);
2854}
2855
2857inline scmatrix operator+(const srmatrix& M1, const scmatrix_slice& M2) {
2858 return spsp_mm_add<srmatrix,scmatrix,scmatrix,complex>(M1,M2.A);
2859}
2860
2862inline scmatrix operator+(const scmatrix& M1, const scmatrix_slice& M2) {
2863 return spsp_mm_add<scmatrix,scmatrix,scmatrix,complex>(M1,M2.A);
2864}
2865
2867inline cmatrix operator+(const scmatrix_slice& M1, const rmatrix& M2) {
2868 return spf_mm_add<scmatrix,rmatrix,cmatrix>(M1.A,M2);
2869}
2870
2872inline cmatrix operator+(const srmatrix_slice& M1, const cmatrix& M2) {
2873 return spf_mm_add<srmatrix,cmatrix,cmatrix>(M1.A,M2);
2874}
2875
2877inline cmatrix operator+(const scmatrix_slice& M1, const cmatrix& M2) {
2878 return spf_mm_add<scmatrix,cmatrix,cmatrix>(M1.A,M2);
2879}
2880
2882inline cmatrix operator+(const cmatrix& M1, const srmatrix_slice& M2) {
2883 return fsp_mm_add<cmatrix,srmatrix,cmatrix>(M1,M2.A);
2884}
2885
2887inline cmatrix operator+(const rmatrix& M1, const scmatrix_slice& M2) {
2888 return fsp_mm_add<rmatrix,scmatrix,cmatrix>(M1,M2.A);
2889}
2890
2892inline cmatrix operator+(const cmatrix& M1, const scmatrix_slice& M2) {
2893 return fsp_mm_add<cmatrix,scmatrix,cmatrix>(M1,M2.A);
2894}
2895
2897inline cmatrix operator+(const scmatrix_slice& M1, const rmatrix_slice& M2) {
2898 return spf_mm_add<scmatrix,rmatrix_slice,cmatrix>(M1.A,M2);
2899}
2900
2902inline cmatrix operator+(const srmatrix_slice& M1, const cmatrix_slice& M2) {
2903 return spf_mm_add<srmatrix,cmatrix_slice,cmatrix>(M1.A,M2);
2904}
2905
2907inline cmatrix operator+(const scmatrix_slice& M1, const cmatrix_slice& M2) {
2908 return spf_mm_add<scmatrix,cmatrix_slice,cmatrix>(M1.A,M2);
2909}
2910
2912inline cmatrix operator+(const cmatrix_slice& M1, const srmatrix_slice& M2) {
2913 return fsp_mm_add<cmatrix_slice,srmatrix,cmatrix>(M1,M2.A);
2914}
2915
2917inline cmatrix operator+(const rmatrix_slice& M1, const scmatrix_slice& M2) {
2918 return fsp_mm_add<rmatrix_slice,scmatrix,cmatrix>(M1,M2.A);
2919}
2920
2922inline cmatrix operator+(const cmatrix_slice& M1, const scmatrix_slice& M2) {
2923 return fsp_mm_add<cmatrix_slice,scmatrix,cmatrix>(M1,M2.A);
2924}
2925
2927inline scmatrix operator-(const scmatrix_slice& M1, const srmatrix_slice& M2) {
2928 return spsp_mm_sub<scmatrix,srmatrix,scmatrix,complex>(M1.A,M2.A);
2929}
2930
2932inline scmatrix operator-(const srmatrix_slice& M1, const scmatrix_slice& M2) {
2933 return spsp_mm_sub<srmatrix,scmatrix,scmatrix,complex>(M1.A,M2.A);
2934}
2935
2937inline scmatrix operator-(const scmatrix_slice& M1, const scmatrix_slice& M2) {
2938 return spsp_mm_sub<scmatrix,scmatrix,scmatrix,complex>(M1.A,M2.A);
2939}
2940
2942inline scmatrix operator-(const scmatrix_slice& M1, const srmatrix& M2) {
2943 return spsp_mm_sub<scmatrix,srmatrix,scmatrix,complex>(M1.A,M2);
2944}
2945
2947inline scmatrix operator-(const srmatrix_slice& M1, const scmatrix& M2) {
2948 return spsp_mm_sub<srmatrix,scmatrix,scmatrix,complex>(M1.A,M2);
2949}
2950
2952inline scmatrix operator-(const scmatrix_slice& M1, const scmatrix& M2) {
2953 return spsp_mm_sub<scmatrix,scmatrix,scmatrix,complex>(M1.A,M2);
2954}
2955
2957inline scmatrix operator-(const scmatrix& M1, const srmatrix_slice& M2) {
2958 return spsp_mm_sub<scmatrix,srmatrix,scmatrix,complex>(M1,M2.A);
2959}
2960
2962inline scmatrix operator-(const srmatrix& M1, const scmatrix_slice& M2) {
2963 return spsp_mm_sub<srmatrix,scmatrix,scmatrix,complex>(M1,M2.A);
2964}
2965
2967inline scmatrix operator-(const scmatrix& M1, const scmatrix_slice& M2) {
2968 return spsp_mm_sub<scmatrix,scmatrix,scmatrix,complex>(M1,M2.A);
2969}
2970
2972inline cmatrix operator-(const scmatrix_slice& M1, const rmatrix& M2) {
2973 return spf_mm_sub<scmatrix,rmatrix,cmatrix>(M1.A,M2);
2974}
2975
2977inline cmatrix operator-(const srmatrix_slice& M1, const cmatrix& M2) {
2978 return spf_mm_sub<srmatrix,cmatrix,cmatrix>(M1.A,M2);
2979}
2980
2982inline cmatrix operator-(const scmatrix_slice& M1, const cmatrix& M2) {
2983 return spf_mm_sub<scmatrix,cmatrix,cmatrix>(M1.A,M2);
2984}
2985
2987inline cmatrix operator-(const cmatrix& M1, const srmatrix_slice& M2) {
2988 return fsp_mm_sub<cmatrix,srmatrix,cmatrix>(M1,M2.A);
2989}
2990
2992inline cmatrix operator-(const rmatrix& M1, const scmatrix_slice& M2) {
2993 return fsp_mm_sub<rmatrix,scmatrix,cmatrix>(M1,M2.A);
2994}
2995
2997inline cmatrix operator-(const cmatrix& M1, const scmatrix_slice& M2) {
2998 return fsp_mm_sub<cmatrix,scmatrix,cmatrix>(M1,M2.A);
2999}
3000
3002inline cmatrix operator-(const scmatrix_slice& M1, const rmatrix_slice& M2) {
3003 return spf_mm_sub<scmatrix,rmatrix_slice,cmatrix>(M1.A,M2);
3004}
3005
3007inline cmatrix operator-(const srmatrix_slice& M1, const cmatrix_slice& M2) {
3008 return spf_mm_sub<srmatrix,cmatrix_slice,cmatrix>(M1.A,M2);
3009}
3010
3012inline cmatrix operator-(const scmatrix_slice& M1, const cmatrix_slice& M2) {
3013 return spf_mm_sub<scmatrix,cmatrix_slice,cmatrix>(M1.A,M2);
3014}
3015
3017inline cmatrix operator-(const cmatrix_slice& M1, const srmatrix_slice& M2) {
3018 return fsp_mm_sub<cmatrix_slice,srmatrix,cmatrix>(M1,M2.A);
3019}
3020
3022inline cmatrix operator-(const rmatrix_slice& M1, const scmatrix_slice& M2) {
3023 return fsp_mm_sub<rmatrix_slice,scmatrix,cmatrix>(M1,M2.A);
3024}
3025
3027inline cmatrix operator-(const cmatrix_slice& M1, const scmatrix_slice& M2) {
3028 return fsp_mm_sub<cmatrix_slice,scmatrix,cmatrix>(M1,M2.A);
3029}
3030
3032 *this = rmatrix(M);
3033 return *this;
3034}
3035
3037 *this = cmatrix(M);
3038 return *this;
3039}
3040
3042 *this += M.A;
3043 return *this;
3044}
3045
3047 *this += M.A;
3048 return *this;
3049}
3050
3052 *this += M.A;
3053 return *this;
3054}
3055
3057 *this += M.A;
3058 return *this;
3059}
3060
3062 *this -= M.A;
3063 return *this;
3064}
3065
3067 *this -= M.A;
3068 return *this;
3069}
3070
3072 *this -= M.A;
3073 return *this;
3074}
3075
3077 *this -= M.A;
3078 return *this;
3079}
3080
3082 *this *= M.A;
3083 return *this;
3084}
3085
3087 *this *= M.A;
3088 return *this;
3089}
3090
3092 *this *= M.A;
3093 return *this;
3094}
3095
3097 *this *= M.A;
3098 return *this;
3099}
3100
3102inline bool operator==(const scmatrix_slice& M1, const srmatrix_slice& M2) {
3103 return spsp_mm_comp(M1.A,M2.A);
3104}
3105
3107inline bool operator==(const srmatrix_slice& M1, const scmatrix_slice& M2) {
3108 return spsp_mm_comp(M1.A,M2.A);
3109}
3110
3112inline bool operator==(const scmatrix_slice& M1, const scmatrix_slice& M2) {
3113 return spsp_mm_comp(M1.A,M2.A);
3114}
3115
3117inline bool operator==(const scmatrix_slice& M1, const srmatrix& M2) {
3118 return spsp_mm_comp(M1.A,M2);
3119}
3120
3122inline bool operator==(const srmatrix_slice& M1, const scmatrix& M2) {
3123 return spsp_mm_comp(M1.A,M2);
3124}
3125
3127inline bool operator==(const scmatrix_slice& M1, const scmatrix& M2) {
3128 return spsp_mm_comp(M1.A,M2);
3129}
3130
3132inline bool operator==(const scmatrix& M1, const srmatrix_slice& M2) {
3133 return spsp_mm_comp(M1,M2.A);
3134}
3135
3137inline bool operator==(const srmatrix& M1, const scmatrix_slice& M2) {
3138 return spsp_mm_comp(M1,M2.A);
3139}
3140
3142inline bool operator==(const scmatrix& M1, const scmatrix_slice& M2) {
3143 return spsp_mm_comp(M1,M2.A);
3144}
3145
3147inline bool operator==(const scmatrix_slice& M1, const rmatrix& M2) {
3148 return spf_mm_comp(M1.A,M2);
3149}
3150
3152inline bool operator==(const srmatrix_slice& M1, const cmatrix& M2) {
3153 return spf_mm_comp(M1.A,M2);
3154}
3155
3157inline bool operator==(const scmatrix_slice& M1, const cmatrix& M2) {
3158 return spf_mm_comp(M1.A,M2);
3159}
3160
3162inline bool operator==(const cmatrix& M1, const srmatrix_slice& M2) {
3163 return fsp_mm_comp(M1,M2.A);
3164}
3165
3167inline bool operator==(const rmatrix& M1, const scmatrix_slice& M2) {
3168 return fsp_mm_comp(M1,M2.A);
3169}
3170
3172inline bool operator==(const cmatrix& M1, const scmatrix_slice& M2) {
3173 return fsp_mm_comp(M1,M2.A);
3174}
3175
3177inline bool operator==(const cmatrix_slice& M1, const srmatrix_slice& M2) {
3178 return fsp_mm_comp(M1,M2.A);
3179}
3180
3182inline bool operator==(const rmatrix_slice& M1, const scmatrix_slice& M2) {
3183 return fsp_mm_comp(M1,M2.A);
3184}
3185
3187inline bool operator==(const cmatrix_slice& M1, const scmatrix_slice& M2) {
3188 return fsp_mm_comp(M1,M2.A);
3189}
3190
3192inline bool operator==(const scmatrix_slice& M1, const rmatrix_slice& M2) {
3193 return spf_mm_comp(M1.A,M2);
3194}
3195
3197inline bool operator==(const srmatrix_slice& M1, const cmatrix_slice& M2) {
3198 return spf_mm_comp(M1.A,M2);
3199}
3200
3202inline bool operator==(const scmatrix_slice& M1, const cmatrix_slice& M2) {
3203 return spf_mm_comp(M1.A,M2);
3204}
3205
3207inline bool operator!=(const scmatrix_slice& M1, const srmatrix_slice& M2) {
3208 return !spsp_mm_comp(M1.A,M2.A);
3209}
3210
3212inline bool operator!=(const srmatrix_slice& M1, const scmatrix_slice& M2) {
3213 return !spsp_mm_comp(M1.A,M2.A);
3214}
3215
3217inline bool operator!=(const scmatrix_slice& M1, const scmatrix_slice& M2) {
3218 return !spsp_mm_comp(M1.A,M2.A);
3219}
3220
3222inline bool operator!=(const scmatrix_slice& M1, const srmatrix& M2) {
3223 return !spsp_mm_comp(M1.A,M2);
3224}
3225
3227inline bool operator!=(const srmatrix_slice& M1, const scmatrix& M2) {
3228 return !spsp_mm_comp(M1.A,M2);
3229}
3230
3232inline bool operator!=(const scmatrix_slice& M1, const scmatrix& M2) {
3233 return !spsp_mm_comp(M1.A,M2);
3234}
3235
3237inline bool operator!=(const scmatrix& M1, const srmatrix_slice& M2) {
3238 return !spsp_mm_comp(M1,M2.A);
3239}
3240
3242inline bool operator!=(const srmatrix& M1, const scmatrix_slice& M2) {
3243 return !spsp_mm_comp(M1,M2.A);
3244}
3245
3247inline bool operator!=(const scmatrix& M1, const scmatrix_slice& M2) {
3248 return !spsp_mm_comp(M1,M2.A);
3249}
3250
3252inline bool operator!=(const scmatrix_slice& M1, const rmatrix& M2) {
3253 return !spf_mm_comp(M1.A,M2);
3254}
3255
3257inline bool operator!=(const srmatrix_slice& M1, const cmatrix& M2) {
3258 return !spf_mm_comp(M1.A,M2);
3259}
3260
3262inline bool operator!=(const scmatrix_slice& M1, const cmatrix& M2) {
3263 return !spf_mm_comp(M1.A,M2);
3264}
3265
3267inline bool operator!=(const cmatrix& M1, const srmatrix_slice& M2) {
3268 return !fsp_mm_comp(M1,M2.A);
3269}
3270
3272inline bool operator!=(const rmatrix& M1, const scmatrix_slice& M2) {
3273 return !fsp_mm_comp(M1,M2.A);
3274}
3275
3277inline bool operator!=(const cmatrix& M1, const scmatrix_slice& M2) {
3278 return !fsp_mm_comp(M1,M2.A);
3279}
3280
3282inline bool operator!=(const cmatrix_slice& M1, const srmatrix_slice& M2) {
3283 return !fsp_mm_comp(M1,M2.A);
3284}
3285
3287inline bool operator!=(const rmatrix_slice& M1, const scmatrix_slice& M2) {
3288 return !fsp_mm_comp(M1,M2.A);
3289}
3290
3292inline bool operator!=(const cmatrix_slice& M1, const scmatrix_slice& M2) {
3293 return !fsp_mm_comp(M1,M2.A);
3294}
3295
3297inline bool operator!=(const scmatrix_slice& M1, const rmatrix_slice& M2) {
3298 return !spf_mm_comp(M1.A,M2);
3299}
3300
3302inline bool operator!=(const srmatrix_slice& M1, const cmatrix_slice& M2) {
3303 return !spf_mm_comp(M1.A,M2);
3304}
3305
3307inline bool operator!=(const scmatrix_slice& M1, const cmatrix_slice& M2) {
3308 return !spf_mm_comp(M1.A,M2);
3309}
3310
3312inline bool operator!(const scmatrix_slice& M) {
3313 return sp_m_not(M.A);
3314}
3315
3317
3322inline std::ostream& operator<<(std::ostream& os, const scmatrix_slice& M) {
3323 return sp_m_output<scmatrix,complex>(os, M.A);
3324}
3325
3327
3332inline std::istream& operator>>(std::istream& is, scmatrix_slice& M) {
3333 scmatrix tmp(M.A.m, M.A.n);
3334 is >> tmp;
3335 M = tmp;
3336 return is;
3337}
3338
3340
3346 private:
3347 scmatrix_slice dat;
3348 bool row;
3349 int index;
3350
3351 scmatrix_subv(scmatrix& A, bool r, int i, int j, int k, int l) : dat(A,i,j,k,l), row(r) {
3352 if(row) index=i; else index=k;
3353 }
3354
3355 scmatrix_subv(const scmatrix& A, bool r, int i, int j, int k, int l) : dat(A,i,j,k,l), row(r) {
3356 if(row) index=i; else index=k;
3357 }
3358
3359
3360 public:
3362
3366 complex& operator[](const int i) {
3367 if(row) {
3368#if(CXSC_INDEX_CHECK)
3369 if(i<dat.A.lb2 || i>dat.A.ub2)
3370 cxscthrow(ELEMENT_NOT_IN_VEC("scmatrix_subv::operator[](int)"));
3371#endif
3372 return dat.element(index,i);
3373 } else {
3374#if(CXSC_INDEX_CHECK)
3375 if(i<dat.A.lb1 || i>dat.A.ub1)
3376 cxscthrow(ELEMENT_NOT_IN_VEC("scmatrix_subv::operator[](int)"));
3377#endif
3378 return dat.element(i,index);
3379 }
3380 }
3381
3383
3386 const complex operator[](const int i) const {
3387 if(row) {
3388#if(CXSC_INDEX_CHECK)
3389 if(i<dat.A.lb2 || i>dat.A.ub2)
3390 cxscthrow(ELEMENT_NOT_IN_VEC("scmatrix_subv::operator[](int)"));
3391#endif
3392 return dat(index,i);
3393 } else {
3394#if(CXSC_INDEX_CHECK)
3395 if(i<dat.A.lb1 || i>dat.A.ub1)
3396 cxscthrow(ELEMENT_NOT_IN_VEC("scmatrix_subv::operator[](int)"));
3397#endif
3398 return dat(i,index);
3399 }
3400 }
3401
3404 return sv_vs_assign(*this,v);
3405 }
3406
3409 return sv_vs_assign(*this,v);
3410 }
3411
3414 return svsp_vv_assign(*this,v);
3415 }
3416
3419 return svsp_vv_assign(*this,v);
3420 }
3421
3424 return svsl_vv_assign(*this,v);
3425 }
3426
3429 return svsl_vv_assign(*this,v);
3430 }
3431
3434 return svf_vv_assign(*this,v);
3435 }
3436
3439 return svf_vv_assign(*this,v);
3440 }
3441
3444 return svf_vv_assign(*this,v);
3445 }
3446
3449 return svf_vv_assign(*this,v);
3450 }
3451
3454 return svsp_vv_assign(*this,srvector(v));
3455 }
3456
3459 return svsp_vv_assign(*this,scvector(v));
3460 }
3461
3462
3464 scmatrix_subv& operator*=(const real&);
3468 scmatrix_subv& operator/=(const real&);
3503
3504 friend scvector operator-(const scmatrix_subv&);
3505 friend std::istream& operator>>(std::istream&, scmatrix_subv&);
3506
3507 friend int Lb(const scmatrix_subv&);
3508 friend int Ub(const scmatrix_subv&);
3509 friend int VecLen(const scmatrix_subv&);
3510 friend srvector Re(const scmatrix_subv&);
3511 friend srvector Im(const scmatrix_subv&);
3512
3513 friend class srvector;
3514 friend class srmatrix;
3515 friend class srmatrix_slice;
3516 friend class scvector;
3517 friend class scmatrix;
3518 friend class scmatrix_slice;
3519 friend class scivector;
3520 friend class scimatrix;
3521 friend class scimatrix_slice;
3522
3523#include "vector_friend_declarations.inl"
3524};
3525
3527inline int Lb(const scmatrix_subv& S) {
3528 if(S.row)
3529 return Lb(S.dat, 2);
3530 else
3531 return Lb(S.dat, 1);
3532}
3533
3535inline int Ub(const scmatrix_subv& S) {
3536 if(S.row)
3537 return Ub(S.dat, 2);
3538 else
3539 return Ub(S.dat, 1);
3540}
3541
3543inline int VecLen(const scmatrix_subv& S) {
3544 return Ub(S)-Lb(S)+1;
3545}
3546
3548inline srvector Re(const scmatrix_subv& S) {
3549 return Re(scvector(S));
3550}
3551
3553inline srvector Im(const scmatrix_subv& S) {
3554 return Im(scvector(S));
3555}
3556
3558inline std::ostream& operator<<(std::ostream& os, const scmatrix_subv& v) {
3559 os << scvector(v);
3560 return os;
3561}
3562
3564inline std::istream& operator>>(std::istream& is, scmatrix_subv& v) {
3565 int n=0;
3566 if(v.row) n=v.dat.A.n; else n=v.dat.A.m;
3567 scvector tmp(n);
3568 is >> tmp;
3569 v = tmp;
3570 return is;
3571}
3572
3573inline scmatrix_subv scmatrix::operator[](const cxscmatrix_column& c) {
3574#if(CXSC_INDEX_CHECK)
3575 if(c.col()<lb2 || c.col()>ub2)
3576 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix::operator[](const cxscmatrix_column&)"));
3577#endif
3578 return scmatrix_subv(*this, false, lb1, ub1, c.col(), c.col());
3579}
3580
3582#if(CXSC_INDEX_CHECK)
3583 if(i<lb1 || i>ub1)
3584 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix::operator[](const int)"));
3585#endif
3586 return scmatrix_subv(*this, true, i, i, lb2, ub2);
3587}
3588
3589inline const scmatrix_subv scmatrix::operator[](const cxscmatrix_column& c) const {
3590#if(CXSC_INDEX_CHECK)
3591 if(c.col()<lb2 || c.col()>ub2)
3592 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix::operator[](const cxscmatrix_column&)"));
3593#endif
3594 return scmatrix_subv(*this, false, lb1, ub1, c.col(), c.col());
3595}
3596
3597inline const scmatrix_subv scmatrix::operator[](const int i) const {
3598#if(CXSC_INDEX_CHECK)
3599 if(i<lb1 || i>ub1)
3600 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix::operator[](const int)"));
3601#endif
3602 return scmatrix_subv(*this, true, i, i, lb2, ub2);
3603}
3604
3606#if(CXSC_INDEX_CHECK)
3607 if(i<A.lb1 || i>A.ub1)
3608 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix_slice::operator[](const int)"));
3609#endif
3610 return scmatrix_subv(*M, true, i, i, A.lb2, A.ub2);
3611}
3612
3613inline scmatrix_subv scmatrix_slice::operator[](const cxscmatrix_column& c) {
3614#if(CXSC_INDEX_CHECK)
3615 if(c.col()<A.lb2 || c.col()>A.ub2)
3616 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix_slice::operator[](const cxscmatrix_column&)"));
3617#endif
3618 return scmatrix_subv(*M, false, A.lb1, A.ub1, c.col(), c.col());
3619}
3620
3621inline const scmatrix_subv scmatrix_slice::operator[](const int i) const {
3622#if(CXSC_INDEX_CHECK)
3623 if(i<A.lb1 || i>A.ub1)
3624 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix_slice::operator[](const int)"));
3625#endif
3626 return scmatrix_subv(*M, true, i, i, A.lb2, A.ub2);
3627}
3628
3629inline const scmatrix_subv scmatrix_slice::operator[](const cxscmatrix_column& c) const {
3630#if(CXSC_INDEX_CHECK)
3631 if(c.col()<A.lb2 || c.col()>A.ub2)
3632 cxscthrow(ROW_OR_COL_NOT_IN_MAT("scmatrix_slice::operator[](const cxscmatrix_column&)"));
3633#endif
3634 return scmatrix_subv(*M, false, A.lb1, A.ub1, c.col(), c.col());
3635}
3636
3638 int nnz = A.dat.A.get_nnz();
3639 p.reserve(nnz);
3640 x.reserve(nnz);
3641
3642 if(A.row) {
3643 lb = A.dat.A.lb2;
3644 ub = A.dat.A.ub2;
3645 n = ub-lb+1;
3646
3647 for(int j=0 ; j<n ; j++) {
3648 for(int k=A.dat.A.p[j] ; k<A.dat.A.p[j+1] ; k++) {
3649 p.push_back(j);
3650 x.push_back(A.dat.A.x[k]);
3651 }
3652 }
3653
3654 } else {
3655 lb = A.dat.A.lb1;
3656 ub = A.dat.A.ub1;
3657 n = ub-lb+1;
3658
3659 for(unsigned int k=0 ; k<A.dat.A.ind.size() ; k++) {
3660 p.push_back(A.dat.A.ind[k]);
3661 x.push_back(A.dat.A.x[k]);
3662 }
3663 }
3664}
3665
3667inline scvector operator-(const scmatrix_subv& v) {
3668 scvector s(v);
3669 return -s;
3670}
3671
3673inline scvector operator/(const scmatrix_subv& v1, const real& v2) {
3674 return scvector(v1) / v2;
3675}
3676
3678inline scvector operator/(const scmatrix_subv& v1, const complex& v2) {
3679 return scvector(v1) / v2;
3680}
3681
3683inline scvector operator/(const srmatrix_subv& v1, const complex& v2) {
3684 return srvector(v1) / v2;
3685}
3686
3688inline scvector operator*(const scmatrix_subv& v1, const real& v2) {
3689 return scvector(v1) * v2;
3690}
3691
3693inline scvector operator*(const scmatrix_subv& v1, const complex& v2) {
3694 return scvector(v1) * v2;
3695}
3696
3698inline scvector operator*(const srmatrix_subv& v1, const complex& v2) {
3699 return srvector(v1) * v2;
3700}
3701
3703inline scvector operator*(const real& v1, const scmatrix_subv& v2) {
3704 return v1 * scvector(v2);
3705}
3706
3708inline scvector operator*(const complex& v1, const scmatrix_subv& v2) {
3709 return v1 * scvector(v2);
3710}
3711
3713inline scvector operator*(const complex& v1, const srmatrix_subv& v2) {
3714 return v1 * srvector(v2);
3715}
3716
3718
3724inline complex operator*(const scmatrix_subv& v1, const srvector& v2) {
3725 return scvector(v1) * v2;
3726}
3727
3729
3735inline complex operator*(const srmatrix_subv& v1, const scvector& v2) {
3736 return srvector(v1) * v2;
3737}
3738
3740
3746inline complex operator*(const scmatrix_subv& v1, const scvector& v2) {
3747 return scvector(v1) * v2;
3748}
3749
3751
3757inline complex operator*(const scmatrix_subv& v1, const srvector_slice& v2) {
3758 return scvector(v1) * v2;
3759}
3760
3762
3768inline complex operator*(const srmatrix_subv& v1, const scvector_slice& v2) {
3769 return srvector(v1) * v2;
3770}
3771
3773
3779inline complex operator*(const scmatrix_subv& v1, const scvector_slice& v2) {
3780 return scvector(v1) * v2;
3781}
3782
3784
3790inline complex operator*(const scmatrix_subv& v1, const rvector& v2) {
3791 return scvector(v1) * v2;
3792}
3793
3795
3801inline complex operator*(const srmatrix_subv& v1, const cvector& v2) {
3802 return srvector(v1) * v2;
3803}
3804
3806
3812inline complex operator*(const scmatrix_subv& v1, const cvector& v2) {
3813 return scvector(v1) * v2;
3814}
3815
3817
3823inline complex operator*(const scmatrix_subv& v1, const rvector_slice& v2) {
3824 return scvector(v1) * v2;
3825}
3826
3828
3834inline complex operator*(const srmatrix_subv& v1, const cvector_slice& v2) {
3835 return srvector(v1) * v2;
3836}
3837
3839
3845inline complex operator*(const scmatrix_subv& v1, const cvector_slice& v2) {
3846 return scvector(v1) * v2;
3847}
3848
3850
3856inline complex operator*(const scvector& v1, const srmatrix_subv& v2) {
3857 return v1 * srvector(v2);
3858}
3859
3861
3867inline complex operator*(const srvector& v1, const scmatrix_subv& v2) {
3868 return v1 * scvector(v2);
3869}
3870
3872
3878inline complex operator*(const scvector& v1, const scmatrix_subv& v2) {
3879 return v1 * scvector(v2);
3880}
3881
3883
3889inline complex operator*(const scvector_slice& v1, const srmatrix_subv& v2) {
3890 return v1 * srvector(v2);
3891}
3892
3894
3900inline complex operator*(const srvector_slice& v1, const scmatrix_subv& v2) {
3901 return v1 * scvector(v2);
3902}
3903
3905
3911inline complex operator*(const scvector_slice& v1, const scmatrix_subv& v2) {
3912 return v1 * scvector(v2);
3913}
3914
3916
3922inline complex operator*(const cvector& v1, const srmatrix_subv& v2) {
3923 return v1 * srvector(v2);
3924}
3925
3927
3933inline complex operator*(const rvector& v1, const scmatrix_subv& v2) {
3934 return v1 * scvector(v2);
3935}
3936
3938
3944inline complex operator*(const cvector& v1, const scmatrix_subv& v2) {
3945 return v1 * scvector(v2);
3946}
3947
3949
3955inline complex operator*(const cvector_slice& v1, const srmatrix_subv& v2) {
3956 return v1 * srvector(v2);
3957}
3958
3960
3966inline complex operator*(const rvector_slice& v1, const scmatrix_subv& v2) {
3967 return v1 * scvector(v2);
3968}
3969
3971
3977inline complex operator*(const cvector_slice& v1, const scmatrix_subv& v2) {
3978 return v1 * scvector(v2);
3979}
3980
3982inline scvector operator+(const scmatrix_subv& v1, const srvector& v2) {
3983 return scvector(v1) + v2;
3984}
3985
3987inline scvector operator+(const srmatrix_subv& v1, const scvector& v2) {
3988 return srvector(v1) + v2;
3989}
3990
3992inline scvector operator+(const scmatrix_subv& v1, const scvector& v2) {
3993 return scvector(v1) + v2;
3994}
3995
3997inline scvector operator+(const scmatrix_subv& v1, const srvector_slice& v2) {
3998 return scvector(v1) + v2;
3999}
4000
4002inline scvector operator+(const srmatrix_subv& v1, const scvector_slice& v2) {
4003 return srvector(v1) + v2;
4004}
4005
4007inline scvector operator+(const scmatrix_subv& v1, const scvector_slice& v2) {
4008 return scvector(v1) + v2;
4009}
4010
4012inline cvector operator+(const scmatrix_subv& v1, const rvector& v2) {
4013 return scvector(v1) + v2;
4014}
4015
4017inline cvector operator+(const srmatrix_subv& v1, const cvector& v2) {
4018 return srvector(v1) + v2;
4019}
4020
4022inline cvector operator+(const scmatrix_subv& v1, const cvector& v2) {
4023 return scvector(v1) + v2;
4024}
4025
4027inline cvector operator+(const scmatrix_subv& v1, const rvector_slice& v2) {
4028 return scvector(v1) + v2;
4029}
4030
4032inline cvector operator+(const srmatrix_subv& v1, const cvector_slice& v2) {
4033 return srvector(v1) + v2;
4034}
4035
4037inline cvector operator+(const scmatrix_subv& v1, const cvector_slice& v2) {
4038 return scvector(v1) + v2;
4039}
4040
4042inline scvector operator+(const scvector& v1, const srmatrix_subv& v2) {
4043 return v1 + srvector(v2);
4044}
4045
4047inline scvector operator+(const srvector& v1, const scmatrix_subv& v2) {
4048 return v1 + scvector(v2);
4049}
4050
4052inline scvector operator+(const scvector& v1, const scmatrix_subv& v2) {
4053 return v1 + scvector(v2);
4054}
4055
4057inline scvector operator+(const scvector_slice& v1, const srmatrix_subv& v2) {
4058 return v1 + srvector(v2);
4059}
4060
4062inline scvector operator+(const srvector_slice& v1, const scmatrix_subv& v2) {
4063 return v1 + scvector(v2);
4064}
4065
4067inline scvector operator+(const scvector_slice& v1, const scmatrix_subv& v2) {
4068 return v1 + scvector(v2);
4069}
4070
4072inline cvector operator+(const cvector& v1, const srmatrix_subv& v2) {
4073 return v1 + srvector(v2);
4074}
4075
4077inline cvector operator+(const rvector& v1, const scmatrix_subv& v2) {
4078 return v1 + scvector(v2);
4079}
4080
4082inline cvector operator+(const cvector& v1, const scmatrix_subv& v2) {
4083 return v1 + scvector(v2);
4084}
4085
4087inline cvector operator+(const cvector_slice& v1, const srmatrix_subv& v2) {
4088 return v1 + srvector(v2);
4089}
4090
4092inline cvector operator+(const rvector_slice& v1, const scmatrix_subv& v2) {
4093 return v1 + scvector(v2);
4094}
4095
4097inline cvector operator+(const cvector_slice& v1, const scmatrix_subv& v2) {
4098 return v1 + scvector(v2);
4099}
4100
4102inline scvector operator-(const scmatrix_subv& v1, const srvector& v2) {
4103 return scvector(v1) - v2;
4104}
4105
4107inline scvector operator-(const srmatrix_subv& v1, const scvector& v2) {
4108 return srvector(v1) - v2;
4109}
4110
4112inline scvector operator-(const scmatrix_subv& v1, const scvector& v2) {
4113 return scvector(v1) - v2;
4114}
4115
4117inline scvector operator-(const scmatrix_subv& v1, const srvector_slice& v2) {
4118 return scvector(v1) - v2;
4119}
4120
4122inline scvector operator-(const srmatrix_subv& v1, const scvector_slice& v2) {
4123 return srvector(v1) - v2;
4124}
4125
4127inline scvector operator-(const scmatrix_subv& v1, const scvector_slice& v2) {
4128 return scvector(v1) - v2;
4129}
4130
4132inline cvector operator-(const scmatrix_subv& v1, const rvector& v2) {
4133 return scvector(v1) - v2;
4134}
4135
4137inline cvector operator-(const srmatrix_subv& v1, const cvector& v2) {
4138 return srvector(v1) - v2;
4139}
4140
4142inline cvector operator-(const scmatrix_subv& v1, const cvector& v2) {
4143 return scvector(v1) - v2;
4144}
4145
4147inline cvector operator-(const scmatrix_subv& v1, const rvector_slice& v2) {
4148 return scvector(v1) - v2;
4149}
4150
4152inline cvector operator-(const srmatrix_subv& v1, const cvector_slice& v2) {
4153 return srvector(v1) - v2;
4154}
4155
4157inline cvector operator-(const scmatrix_subv& v1, const cvector_slice& v2) {
4158 return scvector(v1) - v2;
4159}
4160
4162inline scvector operator-(const scvector& v1, const srmatrix_subv& v2) {
4163 return v1 - srvector(v2);
4164}
4165
4167inline scvector operator-(const srvector& v1, const scmatrix_subv& v2) {
4168 return v1 - scvector(v2);
4169}
4170
4172inline scvector operator-(const scvector& v1, const scmatrix_subv& v2) {
4173 return v1 - scvector(v2);
4174}
4175
4177inline scvector operator-(const scvector_slice& v1, const srmatrix_subv& v2) {
4178 return v1 - srvector(v2);
4179}
4180
4182inline scvector operator-(const srvector_slice& v1, const scmatrix_subv& v2) {
4183 return v1 - scvector(v2);
4184}
4185
4187inline scvector operator-(const scvector_slice& v1, const scmatrix_subv& v2) {
4188 return v1 - scvector(v2);
4189}
4190
4192inline cvector operator-(const cvector& v1, const srmatrix_subv& v2) {
4193 return v1 - srvector(v2);
4194}
4195
4197inline cvector operator-(const rvector& v1, const scmatrix_subv& v2) {
4198 return v1 - scvector(v2);
4199}
4200
4202inline cvector operator-(const cvector& v1, const scmatrix_subv& v2) {
4203 return v1 - scvector(v2);
4204}
4205
4207inline cvector operator-(const cvector_slice& v1, const srmatrix_subv& v2) {
4208 return v1 - srvector(v2);
4209}
4210
4212inline cvector operator-(const rvector_slice& v1, const scmatrix_subv& v2) {
4213 return v1 - scvector(v2);
4214}
4215
4217inline cvector operator-(const cvector_slice& v1, const scmatrix_subv& v2) {
4218 return v1 - scvector(v2);
4219}
4220
4222 *this = *this * v;
4223 return *this;
4224}
4225
4227 *this = *this * v;
4228 return *this;
4229}
4230
4232 *this = *this / v;
4233 return *this;
4234}
4235
4237 *this = *this / v;
4238 return *this;
4239}
4240
4242 *this = *this + v;
4243 return *this;
4244}
4245
4247 *this = *this + v;
4248 return *this;
4249}
4250
4252 *this = *this + v;
4253 return *this;
4254}
4255
4257 *this = *this + v;
4258 return *this;
4259}
4260
4262 *this = *this - v;
4263 return *this;
4264}
4265
4267 *this = *this - v;
4268 return *this;
4269}
4270
4272 *this = *this - v;
4273 return *this;
4274}
4275
4277 *this = *this - v;
4278 return *this;
4279}
4280
4282 *this = *this + v;
4283 return *this;
4284}
4285
4287 *this = *this + v;
4288 return *this;
4289}
4290
4292 *this = *this + v;
4293 return *this;
4294}
4295
4297 *this = *this + v;
4298 return *this;
4299}
4300
4302 *this = *this - v;
4303 return *this;
4304}
4305
4307 *this = *this - v;
4308 return *this;
4309}
4310
4312 *this = *this - v;
4313 return *this;
4314}
4315
4317 *this = *this - v;
4318 return *this;
4319}
4320
4322 *this += rvector(v);
4323 return *this;
4324}
4325
4327 *this += cvector(v);
4328 return *this;
4329}
4330
4332 *this += rvector(v);
4333 return *this;
4334}
4335
4337 *this += cvector(v);
4338 return *this;
4339}
4340
4342 *this += rvector(v);
4343 return *this;
4344}
4345
4347 *this += cvector(v);
4348 return *this;
4349}
4350
4352 *this -= rvector(v);
4353 return *this;
4354}
4355
4357 *this -= cvector(v);
4358 return *this;
4359}
4360
4362 *this = rvector(v);
4363 return *this;
4364}
4365
4367 *this = cvector(v);
4368 return *this;
4369}
4370
4372 *this = rvector(v);
4373 return *this;
4374}
4375
4377 *this = cvector(v);
4378 return *this;
4379}
4380
4382 *this = rvector(v);
4383 return *this;
4384}
4385
4387 *this = cvector(v);
4388 return *this;
4389}
4390
4392inline bool operator==(const scmatrix_subv& v1, const srvector& v2) {
4393 return scvector(v1) == v2;
4394}
4395
4397inline bool operator==(const srmatrix_subv& v1, const scvector& v2) {
4398 return srvector(v1) == v2;
4399}
4400
4402inline bool operator==(const scmatrix_subv& v1, const scvector& v2) {
4403 return scvector(v1) == v2;
4404}
4405
4407inline bool operator==(const scmatrix_subv& v1, const srvector_slice& v2) {
4408 return scvector(v1) == v2;
4409}
4410
4412inline bool operator==(const srmatrix_subv& v1, const scvector_slice& v2) {
4413 return srvector(v1) == v2;
4414}
4415
4417inline bool operator==(const scmatrix_subv& v1, const scvector_slice& v2) {
4418 return scvector(v1) == v2;
4419}
4420
4422inline bool operator==(const scmatrix_subv& v1, const rvector& v2) {
4423 return scvector(v1) == v2;
4424}
4425
4427inline bool operator==(const srmatrix_subv& v1, const cvector& v2) {
4428 return srvector(v1) == v2;
4429}
4430
4432inline bool operator==(const scmatrix_subv& v1, const cvector& v2) {
4433 return scvector(v1) == v2;
4434}
4435
4437inline bool operator==(const scmatrix_subv& v1, const rvector_slice& v2) {
4438 return scvector(v1) == v2;
4439}
4440
4442inline bool operator==(const srmatrix_subv& v1, const cvector_slice& v2) {
4443 return srvector(v1) == v2;
4444}
4445
4447inline bool operator==(const scmatrix_subv& v1, const cvector_slice& v2) {
4448 return scvector(v1) == v2;
4449}
4450
4452inline bool operator==(const scvector& v1, const srmatrix_subv& v2) {
4453 return v1 == srvector(v2);
4454}
4455
4457inline bool operator==(const srvector& v1, const scmatrix_subv& v2) {
4458 return v1 == scvector(v2);
4459}
4460
4462inline bool operator==(const scvector& v1, const scmatrix_subv& v2) {
4463 return v1 == scvector(v2);
4464}
4465
4467inline bool operator==(const scvector_slice& v1, const srmatrix_subv& v2) {
4468 return v1 == srvector(v2);
4469}
4470
4472inline bool operator==(const srvector_slice& v1, const scmatrix_subv& v2) {
4473 return v1 == scvector(v2);
4474}
4475
4477inline bool operator==(const scvector_slice& v1, const scmatrix_subv& v2) {
4478 return v1 == scvector(v2);
4479}
4480
4482inline bool operator==(const cvector& v1, const srmatrix_subv& v2) {
4483 return v1 == srvector(v2);
4484}
4485
4487inline bool operator==(const rvector& v1, const scmatrix_subv& v2) {
4488 return v1 == scvector(v2);
4489}
4490
4492inline bool operator==(const cvector& v1, const scmatrix_subv& v2) {
4493 return v1 == scvector(v2);
4494}
4495
4497inline bool operator==(const cvector_slice& v1, const srmatrix_subv& v2) {
4498 return v1 == srvector(v2);
4499}
4500
4502inline bool operator==(const rvector_slice& v1, const scmatrix_subv& v2) {
4503 return v1 == scvector(v2);
4504}
4505
4507inline bool operator==(const cvector_slice& v1, const scmatrix_subv& v2) {
4508 return v1 == scvector(v2);
4509}
4510
4512inline bool operator!=(const scmatrix_subv& v1, const srvector& v2) {
4513 return scvector(v1) != v2;
4514}
4515
4517inline bool operator!=(const srmatrix_subv& v1, const scvector& v2) {
4518 return srvector(v1) != v2;
4519}
4520
4522inline bool operator!=(const scmatrix_subv& v1, const scvector& v2) {
4523 return scvector(v1) != v2;
4524}
4525
4527inline bool operator!=(const scmatrix_subv& v1, const srvector_slice& v2) {
4528 return scvector(v1) != v2;
4529}
4530
4532inline bool operator!=(const srmatrix_subv& v1, const scvector_slice& v2) {
4533 return srvector(v1) != v2;
4534}
4535
4537inline bool operator!=(const scmatrix_subv& v1, const scvector_slice& v2) {
4538 return scvector(v1) != v2;
4539}
4540
4542inline bool operator!=(const scmatrix_subv& v1, const rvector& v2) {
4543 return scvector(v1) != v2;
4544}
4545
4547inline bool operator!=(const srmatrix_subv& v1, const cvector& v2) {
4548 return srvector(v1) != v2;
4549}
4550
4552inline bool operator!=(const scmatrix_subv& v1, const cvector& v2) {
4553 return scvector(v1) != v2;
4554}
4555
4557inline bool operator!=(const scmatrix_subv& v1, const rvector_slice& v2) {
4558 return scvector(v1) != v2;
4559}
4561
4562inline bool operator!=(const srmatrix_subv& v1, const cvector_slice& v2) {
4563 return srvector(v1) != v2;
4564}
4565
4567inline bool operator!=(const scmatrix_subv& v1, const cvector_slice& v2) {
4568 return scvector(v1) != v2;
4569}
4570
4572inline bool operator!=(const scvector& v1, const srmatrix_subv& v2) {
4573 return v1 != srvector(v2);
4574}
4575
4577inline bool operator!=(const srvector& v1, const scmatrix_subv& v2) {
4578 return v1 != scvector(v2);
4579}
4580
4582inline bool operator!=(const scvector& v1, const scmatrix_subv& v2) {
4583 return v1 != scvector(v2);
4584}
4585
4587inline bool operator!=(const scvector_slice& v1, const srmatrix_subv& v2) {
4588 return v1 != srvector(v2);
4589}
4590
4592inline bool operator!=(const srvector_slice& v1, const scmatrix_subv& v2) {
4593 return v1 != scvector(v2);
4594}
4595
4597inline bool operator!=(const scvector_slice& v1, const scmatrix_subv& v2) {
4598 return v1 != scvector(v2);
4599}
4600
4602inline bool operator!=(const cvector& v1, const srmatrix_subv& v2) {
4603 return v1 != srvector(v2);
4604}
4605
4607inline bool operator!=(const rvector& v1, const scmatrix_subv& v2) {
4608 return v1 != scvector(v2);
4609}
4610
4612inline bool operator!=(const cvector& v1, const scmatrix_subv& v2) {
4613 return v1 != scvector(v2);
4614}
4615
4617inline bool operator!=(const cvector_slice& v1, const srmatrix_subv& v2) {
4618 return v1 != srvector(v2);
4619}
4620
4622inline bool operator!=(const rvector_slice& v1, const scmatrix_subv& v2) {
4623 return v1 != scvector(v2);
4624}
4625
4627inline bool operator!=(const cvector_slice& v1, const scmatrix_subv& v2) {
4628 return v1 != scvector(v2);
4629}
4630
4632inline bool operator!(const scmatrix_subv& x) {
4633 return sv_v_not(x);
4634}
4635
4637
4640inline void accumulate(cdotprecision& dot, const scmatrix_subv& v1, const scmatrix_subv& v2) {
4641 spsp_vv_accu<cdotprecision,scvector,scvector,sparse_cdot>(dot, scvector(v1), scvector(v2));
4642}
4643
4645
4648inline void accumulate(cdotprecision& dot, const scmatrix_subv& v1, const srmatrix_subv& v2) {
4649 spsp_vv_accu<cdotprecision,scvector,srvector,sparse_cdot>(dot, scvector(v1), srvector(v2));
4650}
4651
4653
4656inline void accumulate(cdotprecision& dot, const srmatrix_subv& v1, const scmatrix_subv& v2) {
4657 spsp_vv_accu<cdotprecision,srvector,scvector,sparse_cdot>(dot, srvector(v1), scvector(v2));
4658}
4659
4661
4664inline void accumulate(cdotprecision& dot, const scmatrix_subv& v1, const scvector& v2) {
4665 spsp_vv_accu<cdotprecision,scvector,scvector,sparse_cdot>(dot, scvector(v1), v2);
4666}
4667
4669
4672inline void accumulate(cdotprecision& dot, const scmatrix_subv& v1, const srvector& v2) {
4673 spsp_vv_accu<cdotprecision,scvector,srvector,sparse_cdot>(dot, scvector(v1), v2);
4674}
4675
4677
4680inline void accumulate(cdotprecision& dot, const srmatrix_subv& v1, const scvector& v2) {
4681 spsp_vv_accu<cdotprecision,srvector,scvector,sparse_cdot>(dot, srvector(v1), v2);
4682}
4683
4685
4688inline void accumulate(cdotprecision& dot, const scmatrix_subv& v1, const scvector_slice& v2) {
4689 spsl_vv_accu<cdotprecision,scvector,scvector_slice,sparse_cdot>(dot, scvector(v1), v2);
4690}
4691
4693
4696inline void accumulate(cdotprecision& dot, const scmatrix_subv& v1, const srvector_slice& v2) {
4697 spsl_vv_accu<cdotprecision,scvector,srvector_slice,sparse_cdot>(dot, scvector(v1), v2);
4698}
4699
4701
4704inline void accumulate(cdotprecision& dot, const srmatrix_subv& v1, const scvector_slice& v2) {
4705 spsl_vv_accu<cdotprecision,srvector,scvector_slice,sparse_cdot>(dot, srvector(v1), v2);
4706}
4707
4709
4712inline void accumulate(cdotprecision& dot, const scmatrix_subv& v1, const cvector& v2) {
4713 spf_vv_accu<cdotprecision,scvector,cvector,sparse_cdot>(dot, scvector(v1), v2);
4714}
4715
4717
4720inline void accumulate(cdotprecision& dot, const scmatrix_subv& v1, const rvector& v2) {
4721 spf_vv_accu<cdotprecision,scvector,rvector,sparse_cdot>(dot, scvector(v1), v2);
4722}
4723
4725
4728inline void accumulate(cdotprecision& dot, const srmatrix_subv& v1, const cvector& v2) {
4729 spf_vv_accu<cdotprecision,srvector,cvector,sparse_cdot>(dot, srvector(v1), v2);
4730}
4731
4733
4736inline void accumulate(cdotprecision& dot, const scmatrix_subv& v1, const cvector_slice& v2) {
4737 spf_vv_accu<cdotprecision,scvector,cvector_slice,sparse_cdot>(dot, scvector(v1), v2);
4738}
4739
4741
4744inline void accumulate(cdotprecision& dot, const scmatrix_subv& v1, const rvector_slice& v2) {
4745 spf_vv_accu<cdotprecision,scvector,rvector_slice,sparse_cdot>(dot, scvector(v1), v2);
4746}
4747
4749
4752inline void accumulate(cdotprecision& dot, const srmatrix_subv& v1, const cvector_slice& v2) {
4753 spf_vv_accu<cdotprecision,srvector,cvector_slice,sparse_cdot>(dot, srvector(v1), v2);
4754}
4755
4757
4760inline void accumulate(cdotprecision& dot, const scvector& v1, const scmatrix_subv& v2) {
4761 spsp_vv_accu<cdotprecision,scvector,scvector,sparse_cdot>(dot, v1, scvector(v2));
4762}
4763
4765
4768inline void accumulate(cdotprecision& dot, const scvector& v1, const srmatrix_subv& v2) {
4769 spsp_vv_accu<cdotprecision,scvector,srvector,sparse_cdot>(dot, v1, srvector(v2));
4770}
4771
4773
4776inline void accumulate(cdotprecision& dot, const srvector& v1, const scmatrix_subv& v2) {
4777 spsp_vv_accu<cdotprecision,srvector,scvector,sparse_cdot>(dot, v1, scvector(v2));
4778}
4779
4781
4784inline void accumulate(cdotprecision& dot, const scvector_slice& v1, const scmatrix_subv& v2) {
4785 slsp_vv_accu<cdotprecision,scvector_slice,scvector,sparse_cdot>(dot, v1, scvector(v2));
4786}
4787
4789
4792inline void accumulate(cdotprecision& dot, const scvector_slice& v1, const srmatrix_subv& v2) {
4793 slsp_vv_accu<cdotprecision,scvector_slice,srvector,sparse_cdot>(dot, v1, srvector(v2));
4794}
4795
4797
4800inline void accumulate(cdotprecision& dot, const srvector_slice& v1, const scmatrix_subv& v2) {
4801 slsp_vv_accu<cdotprecision,srvector_slice,scvector,sparse_cdot>(dot, v1, scvector(v2));
4802}
4803
4805
4808inline void accumulate(cdotprecision& dot, const cvector& v1, const scmatrix_subv& v2) {
4809 fsp_vv_accu<cdotprecision,cvector,scvector,sparse_cdot>(dot, v1, scvector(v2));
4810}
4811
4813
4816inline void accumulate(cdotprecision& dot, const cvector& v1, const srmatrix_subv& v2) {
4817 fsp_vv_accu<cdotprecision,cvector,srvector,sparse_cdot>(dot, v1, srvector(v2));
4818}
4819
4821
4824inline void accumulate(cdotprecision& dot, const rvector& v1, const scmatrix_subv& v2) {
4825 fsp_vv_accu<cdotprecision,rvector,scvector,sparse_cdot>(dot, v1, scvector(v2));
4826}
4827
4829
4832inline void accumulate(cdotprecision& dot, const cvector_slice& v1, const scmatrix_subv& v2) {
4833 fsp_vv_accu<cdotprecision,cvector_slice,scvector,sparse_cdot>(dot, v1, scvector(v2));
4834}
4835
4837
4840inline void accumulate(cdotprecision& dot, const cvector_slice& v1, const srmatrix_subv& v2) {
4841 fsp_vv_accu<cdotprecision,cvector_slice,srvector,sparse_cdot>(dot, v1, srvector(v2));
4842}
4843
4845
4848inline void accumulate(cdotprecision& dot, const rvector_slice& v1, const scmatrix_subv& v2) {
4849 fsp_vv_accu<cdotprecision,rvector_slice,scvector,sparse_cdot>(dot, v1, scvector(v2));
4850}
4851
4853
4858inline void accumulate_approx(cdotprecision& dot, const scmatrix_subv& v1, const scmatrix_subv& v2) {
4859 spsp_vv_accuapprox<cdotprecision,scvector,scvector,sparse_cdot>(dot, scvector(v1), scvector(v2));
4860}
4861
4863
4868inline void accumulate_approx(cdotprecision& dot, const scmatrix_subv& v1, const srmatrix_subv& v2) {
4869 spsp_vv_accuapprox<cdotprecision,scvector,srvector,sparse_cdot>(dot, scvector(v1), srvector(v2));
4870}
4871
4873
4878inline void accumulate_approx(cdotprecision& dot, const srmatrix_subv& v1, const scmatrix_subv& v2) {
4879 spsp_vv_accuapprox<cdotprecision,srvector,scvector,sparse_cdot>(dot, srvector(v1), scvector(v2));
4880}
4881
4883
4888inline void accumulate_approx(cdotprecision& dot, const scmatrix_subv& v1, const scvector& v2) {
4889 spsp_vv_accuapprox<cdotprecision,scvector,scvector,sparse_cdot>(dot, scvector(v1), v2);
4890}
4891
4893
4898inline void accumulate_approx(cdotprecision& dot, const scmatrix_subv& v1, const srvector& v2) {
4899 spsp_vv_accuapprox<cdotprecision,scvector,srvector,sparse_cdot>(dot, scvector(v1), v2);
4900}
4901
4903
4908inline void accumulate_approx(cdotprecision& dot, const srmatrix_subv& v1, const scvector& v2) {
4909 spsp_vv_accuapprox<cdotprecision,srvector,scvector,sparse_cdot>(dot, srvector(v1), v2);
4910}
4911
4913
4918inline void accumulate_approx(cdotprecision& dot, const scmatrix_subv& v1, const scvector_slice& v2) {
4919 spsl_vv_accuapprox<cdotprecision,scvector,scvector_slice,sparse_cdot>(dot, scvector(v1), v2);
4920}
4921
4923
4928inline void accumulate_approx(cdotprecision& dot, const scmatrix_subv& v1, const srvector_slice& v2) {
4929 spsl_vv_accuapprox<cdotprecision,scvector,srvector_slice,sparse_cdot>(dot, scvector(v1), v2);
4930}
4931
4933
4938inline void accumulate_approx(cdotprecision& dot, const srmatrix_subv& v1, const scvector_slice& v2) {
4939 spsl_vv_accuapprox<cdotprecision,srvector,scvector_slice,sparse_cdot>(dot, srvector(v1), v2);
4940}
4941
4943
4948inline void accumulate_approx(cdotprecision& dot, const scmatrix_subv& v1, const cvector& v2) {
4949 spf_vv_accuapprox<cdotprecision,scvector,cvector,sparse_cdot>(dot, scvector(v1), v2);
4950}
4951
4953
4958inline void accumulate_approx(cdotprecision& dot, const scmatrix_subv& v1, const rvector& v2) {
4959 spf_vv_accuapprox<cdotprecision,scvector,rvector,sparse_cdot>(dot, scvector(v1), v2);
4960}
4961
4963
4968inline void accumulate_approx(cdotprecision& dot, const srmatrix_subv& v1, const cvector& v2) {
4969 spf_vv_accuapprox<cdotprecision,srvector,cvector,sparse_cdot>(dot, srvector(v1), v2);
4970}
4971
4973
4978inline void accumulate_approx(cdotprecision& dot, const scmatrix_subv& v1, const cvector_slice& v2) {
4979 spf_vv_accuapprox<cdotprecision,scvector,cvector_slice,sparse_cdot>(dot, scvector(v1), v2);
4980}
4981
4983
4988inline void accumulate_approx(cdotprecision& dot, const scmatrix_subv& v1, const rvector_slice& v2) {
4989 spf_vv_accuapprox<cdotprecision,scvector,rvector_slice,sparse_cdot>(dot, scvector(v1), v2);
4990}
4991
4993
4998inline void accumulate_approx(cdotprecision& dot, const srmatrix_subv& v1, const cvector_slice& v2) {
4999 spf_vv_accuapprox<cdotprecision,srvector,cvector_slice,sparse_cdot>(dot, srvector(v1), v2);
5000}
5001
5003
5008inline void accumulate_approx(cdotprecision& dot, const scvector& v1, const scmatrix_subv& v2) {
5009 spsp_vv_accuapprox<cdotprecision,scvector,scvector,sparse_cdot>(dot, v1, scvector(v2));
5010}
5011
5013
5018inline void accumulate_approx(cdotprecision& dot, const scvector& v1, const srmatrix_subv& v2) {
5019 spsp_vv_accuapprox<cdotprecision,scvector,srvector,sparse_cdot>(dot, v1, srvector(v2));
5020}
5021
5023
5028inline void accumulate_approx(cdotprecision& dot, const srvector& v1, const scmatrix_subv& v2) {
5029 spsp_vv_accuapprox<cdotprecision,srvector,scvector,sparse_cdot>(dot, v1, scvector(v2));
5030}
5031
5033
5038inline void accumulate_approx(cdotprecision& dot, const scvector_slice& v1, const scmatrix_subv& v2) {
5039 slsp_vv_accuapprox<cdotprecision,scvector_slice,scvector,sparse_cdot>(dot, v1, scvector(v2));
5040}
5041
5043
5048inline void accumulate_approx(cdotprecision& dot, const scvector_slice& v1, const srmatrix_subv& v2) {
5049 slsp_vv_accuapprox<cdotprecision,scvector_slice,srvector,sparse_cdot>(dot, v1, srvector(v2));
5050}
5051
5053
5058inline void accumulate_approx(cdotprecision& dot, const srvector_slice& v1, const scmatrix_subv& v2) {
5059 slsp_vv_accuapprox<cdotprecision,srvector_slice,scvector,sparse_cdot>(dot, v1, scvector(v2));
5060}
5061
5063
5068inline void accumulate_approx(cdotprecision& dot, const cvector& v1, const scmatrix_subv& v2) {
5069 fsp_vv_accuapprox<cdotprecision,cvector,scvector,sparse_cdot>(dot, v1, scvector(v2));
5070}
5071
5073
5078inline void accumulate_approx(cdotprecision& dot, const cvector& v1, const srmatrix_subv& v2) {
5079 fsp_vv_accuapprox<cdotprecision,cvector,srvector,sparse_cdot>(dot, v1, srvector(v2));
5080}
5081
5083
5088inline void accumulate_approx(cdotprecision& dot, const rvector& v1, const scmatrix_subv& v2) {
5089 fsp_vv_accuapprox<cdotprecision,rvector,scvector,sparse_cdot>(dot, v1, scvector(v2));
5090}
5091
5093
5098inline void accumulate_approx(cdotprecision& dot, const cvector_slice& v1, const scmatrix_subv& v2) {
5099 fsp_vv_accuapprox<cdotprecision,cvector_slice,scvector,sparse_cdot>(dot, v1, scvector(v2));
5100}
5101
5103
5108inline void accumulate_approx(cdotprecision& dot, const cvector_slice& v1, const srmatrix_subv& v2) {
5109 fsp_vv_accuapprox<cdotprecision,cvector_slice,srvector,sparse_cdot>(dot, v1, srvector(v2));
5110}
5111
5113
5118inline void accumulate_approx(cdotprecision& dot, const rvector_slice& v1, const scmatrix_subv& v2) {
5119 fsp_vv_accuapprox<cdotprecision,rvector_slice,scvector,sparse_cdot>(dot, v1, scvector(v2));
5120}
5121
5123
5126inline void accumulate(cidotprecision& dot, const scmatrix_subv& v1, const scmatrix_subv& v2) {
5127 cdotprecision tmp(0.0);
5128 tmp.set_k(dot.get_k());
5129 accumulate(tmp,scvector(v1),scvector(v2));
5130 dot += tmp;
5131}
5132
5134
5137inline void accumulate(cidotprecision& dot, const scmatrix_subv& v1, const srmatrix_subv& v2) {
5138 cdotprecision tmp(0.0);
5139 tmp.set_k(dot.get_k());
5140 accumulate(tmp,scvector(v1),srvector(v2));
5141 dot += tmp;
5142}
5143
5145
5148inline void accumulate(cidotprecision& dot, const srmatrix_subv& v1, const scmatrix_subv& v2) {
5149 cdotprecision tmp(0.0);
5150 tmp.set_k(dot.get_k());
5151 accumulate(tmp,srvector(v1),scvector(v2));
5152 dot += tmp;
5153}
5154
5156
5159inline void accumulate(cidotprecision& dot, const scmatrix_subv& v1, const scvector& v2) {
5160 cdotprecision tmp(0.0);
5161 tmp.set_k(dot.get_k());
5162 accumulate(tmp,scvector(v1),v2);
5163 dot += tmp;
5164}
5165
5167
5170inline void accumulate(cidotprecision& dot, const scmatrix_subv& v1, const srvector& v2) {
5171 cdotprecision tmp(0.0);
5172 tmp.set_k(dot.get_k());
5173 accumulate(tmp,scvector(v1),v2);
5174 dot += tmp;
5175}
5176
5178
5181inline void accumulate(cidotprecision& dot, const srmatrix_subv& v1, const scvector& v2) {
5182 cdotprecision tmp(0.0);
5183 tmp.set_k(dot.get_k());
5184 accumulate(tmp,srvector(v1),v2);
5185 dot += tmp;
5186}
5187
5189
5192inline void accumulate(cidotprecision& dot, const scmatrix_subv& v1, const scvector_slice& v2) {
5193 cdotprecision tmp(0.0);
5194 tmp.set_k(dot.get_k());
5195 accumulate(tmp,scvector(v1),v2);
5196 dot += tmp;
5197}
5198
5200
5203inline void accumulate(cidotprecision& dot, const scmatrix_subv& v1, const srvector_slice& v2) {
5204 cdotprecision tmp(0.0);
5205 tmp.set_k(dot.get_k());
5206 accumulate(tmp,scvector(v1),v2);
5207 dot += tmp;
5208}
5209
5211
5214inline void accumulate(cidotprecision& dot, const srmatrix_subv& v1, const scvector_slice& v2) {
5215 cdotprecision tmp(0.0);
5216 tmp.set_k(dot.get_k());
5217 accumulate(tmp,srvector(v1),v2);
5218 dot += tmp;
5219}
5220
5222
5225inline void accumulate(cidotprecision& dot, const scmatrix_subv& v1, const cvector& v2) {
5226 cdotprecision tmp(0.0);
5227 tmp.set_k(dot.get_k());
5228 accumulate(tmp,scvector(v1),v2);
5229 dot += tmp;
5230}
5231
5233
5236inline void accumulate(cidotprecision& dot, const scmatrix_subv& v1, const rvector& v2) {
5237 cdotprecision tmp(0.0);
5238 tmp.set_k(dot.get_k());
5239 accumulate(tmp,scvector(v1),v2);
5240 dot += tmp;
5241}
5242
5244
5247inline void accumulate(cidotprecision& dot, const srmatrix_subv& v1, const cvector& v2) {
5248 cdotprecision tmp(0.0);
5249 tmp.set_k(dot.get_k());
5250 accumulate(tmp,srvector(v1),v2);
5251 dot += tmp;
5252}
5253
5255
5258inline void accumulate(cidotprecision& dot, const scmatrix_subv& v1, const cvector_slice& v2) {
5259 cdotprecision tmp(0.0);
5260 tmp.set_k(dot.get_k());
5261 accumulate(tmp,scvector(v1),v2);
5262 dot += tmp;
5263}
5264
5266
5269inline void accumulate(cidotprecision& dot, const scmatrix_subv& v1, const rvector_slice& v2) {
5270 cdotprecision tmp(0.0);
5271 tmp.set_k(dot.get_k());
5272 accumulate(tmp,scvector(v1),v2);
5273 dot += tmp;
5274}
5275
5277
5280inline void accumulate(cidotprecision& dot, const srmatrix_subv& v1, const cvector_slice& v2) {
5281 cdotprecision tmp(0.0);
5282 tmp.set_k(dot.get_k());
5283 accumulate(tmp,srvector(v1),v2);
5284 dot += tmp;
5285}
5286
5288
5291inline void accumulate(cidotprecision& dot, const scvector& v1, const scmatrix_subv& v2) {
5292 cdotprecision tmp(0.0);
5293 tmp.set_k(dot.get_k());
5294 accumulate(tmp,v1,scvector(v2));
5295 dot += tmp;
5296}
5297
5299
5302inline void accumulate(cidotprecision& dot, const scvector& v1, const srmatrix_subv& v2) {
5303 cdotprecision tmp(0.0);
5304 tmp.set_k(dot.get_k());
5305 accumulate(tmp,v1,srvector(v2));
5306 dot += tmp;
5307}
5308
5310
5313inline void accumulate(cidotprecision& dot, const srvector& v1, const scmatrix_subv& v2) {
5314 cdotprecision tmp(0.0);
5315 tmp.set_k(dot.get_k());
5316 accumulate(tmp,v1,scvector(v2));
5317 dot += tmp;
5318}
5319
5321
5324inline void accumulate(cidotprecision& dot, const scvector_slice& v1, const scmatrix_subv& v2) {
5325 cdotprecision tmp(0.0);
5326 tmp.set_k(dot.get_k());
5327 accumulate(tmp,v1,scvector(v2));
5328 dot += tmp;
5329}
5330
5332
5335inline void accumulate(cidotprecision& dot, const scvector_slice& v1, const srmatrix_subv& v2) {
5336 cdotprecision tmp(0.0);
5337 tmp.set_k(dot.get_k());
5338 accumulate(tmp,v1,srvector(v2));
5339 dot += tmp;
5340}
5341
5343
5346inline void accumulate(cidotprecision& dot, const srvector_slice& v1, const scmatrix_subv& v2) {
5347 cdotprecision tmp(0.0);
5348 tmp.set_k(dot.get_k());
5349 accumulate(tmp,v1,scvector(v2));
5350 dot += tmp;
5351}
5352
5354
5357inline void accumulate(cidotprecision& dot, const cvector& v1, const scmatrix_subv& v2) {
5358 cdotprecision tmp(0.0);
5359 tmp.set_k(dot.get_k());
5360 accumulate(tmp,v1,scvector(v2));
5361 dot += tmp;
5362}
5363
5365
5368inline void accumulate(cidotprecision& dot, const cvector& v1, const srmatrix_subv& v2) {
5369 cdotprecision tmp(0.0);
5370 tmp.set_k(dot.get_k());
5371 accumulate(tmp,v1,srvector(v2));
5372 dot += tmp;
5373}
5374
5376
5379inline void accumulate(cidotprecision& dot, const rvector& v1, const scmatrix_subv& v2) {
5380 cdotprecision tmp(0.0);
5381 tmp.set_k(dot.get_k());
5382 accumulate(tmp,v1,scvector(v2));
5383 dot += tmp;
5384}
5385
5387
5390inline void accumulate(cidotprecision& dot, const cvector_slice& v1, const scmatrix_subv& v2) {
5391 cdotprecision tmp(0.0);
5392 tmp.set_k(dot.get_k());
5393 accumulate(tmp,v1,scvector(v2));
5394 dot += tmp;
5395}
5396
5398
5401inline void accumulate(cidotprecision& dot, const cvector_slice& v1, const srmatrix_subv& v2) {
5402 cdotprecision tmp(0.0);
5403 tmp.set_k(dot.get_k());
5404 accumulate(tmp,v1,srvector(v2));
5405 dot += tmp;
5406}
5407
5409
5412inline void accumulate(cidotprecision& dot, const rvector_slice& v1, const scmatrix_subv& v2) {
5413 cdotprecision tmp(0.0);
5414 tmp.set_k(dot.get_k());
5415 accumulate(tmp,v1,scvector(v2));
5416 dot += tmp;
5417}
5418
5419} //namespace cxsc;
5420
5421#include "sparsematrix.inl"
5422
5423#endif
5424
The Data Type cdotprecision.
Definition cdot.hpp:61
void set_k(unsigned int i)
Set precision for computation of dot products.
Definition cdot.hpp:93
The Data Type cidotprecision.
Definition cidot.hpp:58
int get_k() const
Get currently set precision for computation of dot products.
Definition cidot.hpp:89
The Data Type cimatrix.
Definition cimatrix.hpp:908
The Data Type cmatrix_slice.
Definition cmatrix.hpp:1203
cmatrix_slice & operator=(const cmatrix &m) noexcept
Implementation of standard assigning operator.
Definition cmatrix.inl:450
cmatrix_slice & operator-=(const srmatrix &m)
Implementation of addition and assignment operator.
cmatrix_slice & operator*=(const srmatrix &m)
Implementation of addition and assignment operator.
cmatrix_slice & operator+=(const srmatrix &m)
Implementation of addition and assignment operator.
The Data Type cmatrix_subv.
Definition cmatrix.hpp:54
cmatrix_subv & operator=(const scmatrix_subv &rv)
Implementation of standard assigning operator.
cmatrix_subv & operator-=(const scmatrix_subv &rv)
Implementation of standard assigning operator.
cmatrix_subv & operator+=(const scmatrix_subv &rv)
Implementation of standard assigning operator.
The Data Type cmatrix.
Definition cmatrix.hpp:514
cmatrix() noexcept
Constructor of class cmatrix.
Definition cmatrix.inl:31
cmatrix & operator+=(const srmatrix &m)
Implementation of addition and assignment operator.
cmatrix & operator*=(const srmatrix &m)
Implementation of addition and assignment operator.
cmatrix & operator-=(const srmatrix &m)
Implementation of addition and assignment operator.
cmatrix & operator=(const complex &r) noexcept
Implementation of standard assigning operator.
Definition cmatrix.inl:438
The Scalar Type complex.
Definition complex.hpp:50
The Data Type cvector_slice.
Definition cvector.hpp:845
The Data Type cvector.
Definition cvector.hpp:58
The Data Type intmatrix.
The Data Type intvector.
Definition intvector.hpp:52
The Scalar Type real.
Definition real.hpp:114
The Data Type rmatrix_slice.
Definition rmatrix.hpp:1443
The Data Type rmatrix.
Definition rmatrix.hpp:471
The Data Type rvector_slice.
Definition rvector.hpp:1064
The Data Type rvector.
Definition rvector.hpp:58
A slice of a sparse complex interval matrix.
Represents a row or column vector of a sparse matrix.
A sparse complex interval matrix.
Definition scimatrix.hpp:71
A sparse complex interval vector.
Definition scivector.hpp:62
A slice of a sparse complex matrix.
scmatrix_slice & operator*=(const cmatrix_slice &M)
Assigns the product of the sparse slice and M to the slice.
scmatrix_slice & operator-=(const cmatrix_slice &M)
Assigns the element wise difference of the sparse slice and M to the slice.
friend int RowLen(const scmatrix_slice &)
Returns the number columns of the matrix slice.
friend int Ub(const scmatrix_slice &, const int)
Returns the upper index bound of the rows (if i==ROW) or columns (if i==COL) of the slice.
scmatrix_slice & operator-=(const scmatrix_slice &M)
Assigns the element wise difference of the sparse slice and M to the slice.
scmatrix_slice & operator*=(const real &r)
Assigns the component wise product of the sparse slice and r to the slice.
scmatrix_slice & operator+=(const cmatrix_slice &M)
Assigns the element wise sum of the sparse slice and M to the slice.
scmatrix_slice & operator-=(const rmatrix &M)
Assigns the element wise difference of the sparse slice and M to the slice.
friend srmatrix Re(const scmatrix_slice &)
Return the real part of the slice.
scmatrix_slice & operator=(const real &C)
Assing C to all elements of the slice
scmatrix_slice & operator+=(const scmatrix &M)
Assigns the element wise sum of the sparse slice and M to the slice.
scmatrix_slice & operator+=(const srmatrix_slice &M)
Assigns the element wise sum of the sparse slice and M to the slice.
scmatrix_slice & operator=(const cmatrix_slice &C)
Assing C to the slice.
scmatrix_subv operator[](const int)
Returns a row of the matrix.
scmatrix_slice & operator/=(const real &r)
Assigns the component wise division of the sparse slice and M to the slice.
scmatrix_slice & operator-=(const srmatrix_slice &M)
Assigns the element wise difference of the sparse slice and M to the slice.
scmatrix_slice & operator=(const rmatrix &C)
Assing C to the slice.
scmatrix_slice & operator/=(const complex &r)
Assigns the component wise division of the sparse slice and M to the slice.
friend std::ostream & operator<<(std::ostream &, const scmatrix_slice &)
Standard output operator for sparse matrix slice.
scmatrix_slice & operator*=(const scmatrix &M)
Assigns the product of the sparse slice and M to the slice.
scmatrix_slice & operator+=(const scmatrix_slice &M)
Assigns the element wise sum of the sparse slice and M to the slice.
scmatrix_slice & operator*=(const srmatrix_slice &M)
Assigns the product of the sparse slice and M to the slice.
scmatrix_slice & operator+=(const rmatrix_slice &M)
Assigns the element wise sum of the sparse slice and M to the slice.
scmatrix_slice & operator=(const cmatrix &C)
Assing C to the slice.
friend int Lb(const scmatrix_slice &, const int)
Returns the lower index bound of the rows (if i==ROW) or columns (if i==COL) of the slice.
scmatrix_slice & operator*=(const rmatrix &M)
Assigns the product of the sparse slice and M to the slice.
scmatrix_slice & operator+=(const cmatrix &M)
Assigns the element wise sum of the sparse slice and M to the slice.
scmatrix_slice & operator=(const complex &C)
Assing C to all elements of the slice.
scmatrix_slice & operator=(const scmatrix &C)
Assing C to the slice.
scmatrix_slice & operator-=(const rmatrix_slice &M)
Assigns the element wise difference of the sparse slice and M to the slice.
scmatrix_slice & operator=(const srmatrix_slice &C)
Assing C to the slice.
friend srmatrix Im(const scmatrix_slice &)
Returns the imaginary part of the slice.
scmatrix_slice & operator*=(const srmatrix &M)
Assigns the product of the sparse slice and M to the slice.
complex & element(const int i, const int j)
Returns a reference to the element (i,j) of the matrix.
scmatrix_slice & operator=(const srmatrix &C)
Assing C to the slice.
scmatrix_slice & operator=(const scmatrix_slice &C)
Assing C to the slice.
scmatrix_slice & operator*=(const rmatrix_slice &M)
Assigns the product of the sparse slice and M to the slice.
scmatrix_slice & operator*=(const cmatrix &M)
Assigns the product of the sparse slice and M to the slice.
scmatrix_slice & operator+=(const srmatrix &M)
Assigns the element wise sum of the sparse slice and M to the slice.
const complex operator()(const int i, const int j) const
Returns a copy of the element (i,j) of the matrix.
scmatrix_slice & operator-=(const scmatrix &M)
Assigns the element wise difference of the sparse slice and M to the slice.
scmatrix_slice & operator+=(const rmatrix &M)
Assigns the element wise sum of the sparse slice and M to the slice.
scmatrix_slice & operator-=(const cmatrix &M)
Assigns the element wise difference of the sparse slice and M to the slice.
scmatrix_slice & operator*=(const scmatrix_slice &M)
Assigns the product of the sparse slice and M to the slice.
scmatrix_slice & operator=(const rmatrix_slice &C)
Assing C to the slice.
friend int ColLen(const scmatrix_slice &)
Returns the number of rows of the matrix slice.
scmatrix_slice & operator-=(const srmatrix &M)
Assigns the element wise difference of the sparse slice and M to the slice.
scmatrix_slice & operator*=(const complex &r)
Assigns the component wise product of the sparse slice and r to the slice.
Represents a row or column vector of a sparse matrix.
scmatrix_subv & operator/=(const real &)
Assign the componentwise division of the subvector with a scalar to the subvector
scmatrix_subv & operator-=(const srvector &)
Assign the difference of the subvector with a vector to the subvector
complex & operator[](const int i)
Returns a reference to the i-th element of the subvector.
scmatrix_subv & operator=(const scvector_slice &v)
Assigns a vector to a subvector.
scmatrix_subv & operator*=(const real &)
Assign the componentwise product of the subvector with a scalar to the subvector.
friend int Lb(const scmatrix_subv &)
Returns the lower index bound of the subvector.
friend scvector operator-(const scmatrix_subv &)
Unary negation operator.
scmatrix_subv & operator=(const srvector &v)
Assigns a vector to a subvector.
scmatrix_subv & operator+=(const srvector &)
Assign the sum of the subvector with a vector to the subvector
scmatrix_subv & operator=(const rvector &v)
Assigns a vector to a subvector.
scmatrix_subv & operator=(const cvector &v)
Assigns a vector to a subvector.
friend std::istream & operator>>(std::istream &, scmatrix_subv &)
Standard input operator for subvectors.
scmatrix_subv & operator=(const srmatrix_subv &v)
Assigns a vector to a subvector.
friend int Ub(const scmatrix_subv &)
Returns the upper index bound of the subvector.
const complex operator[](const int i) const
Returns a copy of the i-th element of the subvector.
scmatrix_subv & operator=(const srvector_slice &v)
Assigns a vector to a subvector.
friend int VecLen(const scmatrix_subv &)
Returns the length of the subvector.
scmatrix_subv & operator=(const real &v)
Assigns v to all elements of the subvector.
scmatrix_subv & operator=(const complex &v)
Assigns v to all elements of the subvector.
scmatrix_subv & operator=(const rvector_slice &v)
Assigns a vector to a subvector.
friend srvector Re(const scmatrix_subv &)
Returns the real part of the subvector.
friend srvector Im(const scmatrix_subv &)
Returns the imaginary part of the subvector.
scmatrix_subv & operator=(const cvector_slice &v)
Assigns a vector to a subvector.
scmatrix_subv & operator=(const scvector &v)
Assigns a vector to a subvector.
scmatrix_subv & operator=(const scmatrix_subv &v)
Assigns a vector to a subvector.
A sparse complex matrix.
Definition scmatrix.hpp:69
friend srmatrix CompMat(const scmatrix &)
Returns Ostrowskis comparison matrix for A.
int get_nnz() const
Returns the number of non-zero entries (including explicitly stored zeros).
Definition scmatrix.hpp:684
scmatrix operator()(const intmatrix &P, const intmatrix &Q)
Performs row and column permutations using the two permutation matrices P and Q. Faster than explicit...
Definition scmatrix.hpp:666
friend int ColLen(const scmatrix &)
Returns the number of rows of the matrix.
Definition scmatrix.hpp:984
scmatrix(const srmatrix &A)
Creates a sparse complex matrix out of a sparse real matrix.
Definition scmatrix.hpp:359
scmatrix & operator*=(const cmatrix_slice &B)
Multiply the sparse matrix by B and assign the result to it.
Definition scmatrix.hpp:764
scmatrix & operator=(const srmatrix &A)
Assigns a sparse real matrix to the sparse complex matrix.
Definition scmatrix.hpp:527
real density() const
Returns the density (the number of non-zeros divided by the number of elements) of the matrix.
Definition scmatrix.hpp:679
const std::vector< int > & row_indices() const
Returns a constant reference to the vector containing the row indices (the array)
Definition scmatrix.hpp:102
friend srmatrix Re(const scmatrix &)
Returns the real part of the sparse matrix A.
friend void SetUb(scmatrix &, const int, const int)
Sets the upper index bound of the rows (i==ROW) or columns (i==COL) to j.
Definition scmatrix.hpp:941
friend scmatrix mid(const scimatrix &)
Returns the componentwise midpoint of the matrix A.
const std::vector< complex > & values() const
Returns a constant reference to the vector containing the stored values (the array)
Definition scmatrix.hpp:107
scmatrix & operator-=(const cmatrix_slice &B)
Subtract B from the sparse matrix and assign the result to it.
Definition scmatrix.hpp:734
scmatrix & operator-=(const cmatrix &B)
Subtract B from the sparse matrix and assign the result to it.
Definition scmatrix.hpp:724
scmatrix(const int m, const int n, const int nnz, const intvector &rows, const intvector &cols, const cvector &values, const enum STORAGE_TYPE t=triplet)
Creates a sparse matrix out of three vectors (arrays) forming a matrix stored in triplet,...
Definition scmatrix.hpp:143
scmatrix & operator*=(const rmatrix &B)
Multiply the sparse matrix by B and assign the result to it.
Definition scmatrix.hpp:754
scmatrix & operator*=(const complex &r)
Multiply all elements of the sparse matrix by r and assign the result to it.
Definition scmatrix.hpp:784
scmatrix & operator+=(const scmatrix &B)
Add B to the sparse matrix and assign the result to it.
Definition scmatrix.hpp:714
scmatrix_subv operator[](const cxscmatrix_column &)
Returns a column of the matrix as a sparse subvector object.
friend scmatrix Sup(const scimatrix &)
Returns the Supremum of the matrix A.
scmatrix & operator+=(const cmatrix &B)
Add B to the sparse matrix and assign the result to it.
Definition scmatrix.hpp:694
scmatrix & operator=(const rmatrix &A)
Assigns a dense matrix to the sparse matrix. Only the non zero entries of the dense matrix slice are ...
Definition scmatrix.hpp:507
scmatrix & operator/=(const real &r)
Divide all elements of the sparse matrix by r and assign the result to it.
Definition scmatrix.hpp:789
scmatrix & operator*=(const srmatrix &B)
Multiply the sparse matrix by B and assign the result to it.
Definition scmatrix.hpp:769
friend scmatrix diam(const scimatrix &)
Returns the componentwise diameter of the matrix A.
friend int Lb(const scmatrix &, int)
Returns the lower index bound for the rows or columns of A.
Definition scmatrix.hpp:956
scmatrix & operator-=(const rmatrix_slice &B)
Subtract B from the sparse matrix and assign the result to it.
Definition scmatrix.hpp:729
scmatrix & operator+=(const srmatrix &B)
Add B to the sparse matrix and assign the result to it.
Definition scmatrix.hpp:709
void dropzeros()
Drops explicitly stored zeros from the data structure.
Definition scmatrix.hpp:473
scmatrix operator()(const intvector &pervec, const intvector &q)
Performs a row and column permutation using two permutation vectors.
Definition scmatrix.hpp:615
scmatrix & operator/=(const complex &r)
Divide all elements of the sparse matrix by r and assign the result to it.
Definition scmatrix.hpp:794
scmatrix & operator-=(const rmatrix &B)
Subtract B from the sparse matrix and assign the result to it.
Definition scmatrix.hpp:719
scmatrix & operator=(const real &A)
Assigns a real value to all elements of the matrix (resulting in a dense matrix!)
Definition scmatrix.hpp:497
scmatrix & operator-=(const scmatrix &B)
Subtract B from the sparse matrix and assign the result to it.
Definition scmatrix.hpp:744
scmatrix & operator=(const cmatrix_slice &A)
Assigns a dense matrix slice to the sparse matrix. Only the non zero entries of the dense matrix slic...
Definition scmatrix.hpp:522
void full(cmatrix &A) const
Creates a full matrix out of the sparse matrix and stores it in A. This should normally be done using...
Definition scmatrix.hpp:458
friend srmatrix abs(const scmatrix &)
Returns the componentwise absolute value of the sparse matrix A.
scmatrix & operator+=(const cmatrix_slice &B)
Add B to the sparse matrix and assign the result to it.
Definition scmatrix.hpp:704
scmatrix & operator*=(const real &r)
Multiply all elements of the sparse matrix by r and assign the result to it.
Definition scmatrix.hpp:779
std::vector< complex > & values()
Returns a reference to the vector containing the stored values (the array)
Definition scmatrix.hpp:92
complex & element(int i, int j)
Returns a reference to the element (i,j) of the matrix.
Definition scmatrix.hpp:579
scmatrix operator()(const intmatrix &P)
Performs a row permutation using the permutation matrix P. Faster than explicitly computing the produ...
Definition scmatrix.hpp:673
friend int RowLen(const scmatrix &)
Returns the number of columns of the matrix.
Definition scmatrix.hpp:979
scmatrix(const int m, const int n, const int nnz, const int *rows, const int *cols, const complex *values, const enum STORAGE_TYPE t=triplet)
Creates a sparse matrix out of three arrays forming a matrix stored in triplet, compressed row or com...
Definition scmatrix.hpp:254
scmatrix & operator=(const cmatrix &A)
Assigns a dense matrix to the sparse matrix. Only the non zero entries of the dense matrix slice are ...
Definition scmatrix.hpp:512
std::vector< int > & row_indices()
Returns a reference to the vector containing the row indices (the array)
Definition scmatrix.hpp:87
scmatrix()
Standard constructor, creates an empty matrix of dimension 0x0.
Definition scmatrix.hpp:112
scmatrix(const cmatrix &A)
Creates a sparse matrix out of a dense matrix A. Only the non zero elements of A are stored explicitl...
Definition scmatrix.hpp:390
friend scmatrix transp(const scmatrix &)
Returns the transpose of A.
Definition scmatrix.hpp:886
friend std::istream & operator>>(std::istream &, scmatrix_slice &)
Standard input operator for sparse matrix slice.
const std::vector< int > & column_pointers() const
Returns a constant reference to the vector containing the column pointers (the array)
Definition scmatrix.hpp:97
scmatrix operator()(const intvector &pervec)
Performs a row permutation using a permutation vector.
Definition scmatrix.hpp:642
scmatrix & operator=(const complex &A)
Assigns a complex value to all elements of the matrix (resulting in a dense matrix!...
Definition scmatrix.hpp:502
scmatrix & operator+=(const rmatrix &B)
Add B to the sparse matrix and assign the result to it.
Definition scmatrix.hpp:689
const complex operator()(int i, int j) const
Returns a copy of the element in row i and column j.
Definition scmatrix.hpp:558
scmatrix(const rmatrix &A)
Creates a sparse matrix out of a dense matrix A. Only the non zero elements of A are stored explicitl...
Definition scmatrix.hpp:367
friend srmatrix Im(const scmatrix &)
Returns the imaginary part of the sparse matrix A.
scmatrix(const int ms, const int ns, const cmatrix &A)
Constructor for banded matrices.
Definition scmatrix.hpp:416
scmatrix & operator*=(const rmatrix_slice &B)
Multiply the sparse matrix by B and assign the result to it.
Definition scmatrix.hpp:759
friend void SetLb(scmatrix &, const int, const int)
Sets the lower index bound of the rows (i==ROW) or columns (i==COL) to j.
Definition scmatrix.hpp:926
scmatrix(const int r, const int c)
Creates an empty matrix with r rows and c columns, pre-reserving space for 2*(r+c) elements.
Definition scmatrix.hpp:119
scmatrix & operator*=(const scmatrix &B)
Multiply the sparse matrix by B and assign the result to it.
Definition scmatrix.hpp:774
scmatrix(const int r, const int c, const int e)
Creates an empty matrix with r rows and c columns, pre-reserving space for e elements.
Definition scmatrix.hpp:128
std::vector< int > & column_pointers()
Returns a reference to the vector containing the column pointers (the array)
Definition scmatrix.hpp:82
scmatrix & operator+=(const rmatrix_slice &B)
Add B to the sparse matrix and assign the result to it.
Definition scmatrix.hpp:699
scmatrix & operator-=(const srmatrix &B)
Subtract B from the sparse matrix and assign the result to it.
Definition scmatrix.hpp:739
scmatrix & operator=(const rmatrix_slice &A)
Assigns a dense matrix slice to the sparse matrix. Only the non zero entries of the dense matrix slic...
Definition scmatrix.hpp:517
friend int Ub(const scmatrix &, int)
Returns the upper index bound for the rows or columns of A.
Definition scmatrix.hpp:969
scmatrix & operator*=(const cmatrix &B)
Multiply the sparse matrix by B and assign the result to it.
Definition scmatrix.hpp:749
friend scmatrix Id(const scmatrix &)
Return a sparse unity matrix of the same dimension as A.
Definition scmatrix.hpp:863
friend scmatrix Inf(const scimatrix &)
Returns the Infimum of the matrix A.
Helper class for slices of sparse vectors.
A sparse complex vector.
Definition scvector.hpp:58
scvector()
Default constructor, creates an empty vector of size 0
Definition scvector.hpp:68
A slice of a sparse real matrix.
Represents a row or column vector of a sparse matrix.
A sparse real matrix.
Definition srmatrix.hpp:77
void dropzeros()
Drops explicitly stored zeros from the data structure.
Definition srmatrix.hpp:449
int get_nnz() const
Returns the number of non-zero entries (including explicitly stored zeros).
Definition srmatrix.hpp:630
Helper class for slices of sparse vectors.
Definition srvector.hpp:868
A sparse real vector.
Definition srvector.hpp:58
The namespace cxsc, providing all functionality of the class library C-XSC.
Definition cdot.cpp:29
int VecLen(const scimatrix_subv &S)
Returns the length of the subvector.
void accumulate_approx(cdotprecision &dp, const cmatrix_subv &rv1, const cmatrix_subv &rv2)
The accurate scalar product of the last two arguments added to the value of the first argument (witho...
Definition cmatrix.cpp:99
int ColLen(const cimatrix &)
Returns the column dimension.
rmatrix CompMat(const cimatrix &A)
Returns Ostrowski's comparison matrix.
Definition cimatrix.cpp:45
civector operator/(const cimatrix_subv &rv, const cinterval &s) noexcept
Implementation of division operation.
Definition cimatrix.inl:730
cimatrix transp(const cimatrix &A)
Returns the transposed matrix.
Definition cimatrix.cpp:74
int Ub(const cimatrix &rm, const int &i) noexcept
Returns the upper bound index.
cimatrix & SetLb(cimatrix &m, const int &i, const int &j) noexcept
Sets the lower bound index.
int RowLen(const cimatrix &)
Returns the row dimension.
STORAGE_TYPE
Enumeration depicting the storage type of a sparse matrix (Triplet storage, Compressed column storage...
Definition srmatrix.hpp:42
cimatrix Id(const cimatrix &A)
Returns the Identity matrix.
Definition cimatrix.cpp:61
ivector abs(const cimatrix_subv &mv) noexcept
Returns the absolute value of the matrix.
Definition cimatrix.inl:737
civector operator*(const cimatrix_subv &rv, const cinterval &s) noexcept
Implementation of multiplication operation.
Definition cimatrix.inl:731
void Resize(cimatrix &A) noexcept
Resizes the matrix.
cimatrix & SetUb(cimatrix &m, const int &i, const int &j) noexcept
Sets the upper bound index.
int Lb(const cimatrix &rm, const int &i) noexcept
Returns the lower bound index.