fmatvec 0.0.0

general_matrix.h
1 /* Copyright (C) 2003-2005 Martin Förg
2
3 * This library is free software; you can redistribute it and/or
4 * modify it under the terms of the GNU Lesser General Public
5 * License as published by the Free Software Foundation; either
6 * version 2.1 of the License, or (at your option) any later version.
7 *
8 * This library is distributed in the hope that it will be useful,
9 * but WITHOUT ANY WARRANTY; without even the implied warranty of
10 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 * Lesser General Public License for more details.
12 *
13 * You should have received a copy of the GNU Lesser General Public
14 * License along with this library; if not, write to the Free Software
15 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
16 *
17 * Contact:
18 * martin.o.foerg@googlemail.com
19 *
20 */
21
22 #ifndef general_matrix_h
23 #define general_matrix_h
24
25 #include "_memory.h"
26 #include "types.h"
27 #include "indices.h"
28 #include "matrix.h"
29 #include <cstdlib>
30 #include <stdexcept>
31
32 namespace fmatvec {
33
42 template <class AT> class Matrix<General,Ref,Ref,AT> {
43
44 public:
45 static constexpr bool isVector {false};
46
47 using value_type = AT;
48 using shape_type = General;
49
51
52 friend class Vector<Ref,AT>;
53 friend class RowVector<Ref,AT>;
54 friend class SquareMatrix<Ref,AT>;
55 friend class Matrix<Symmetric,Ref,Ref,AT>;
56
57 protected:
58
59 template<bool Const>
60 class Iterator {
61 using Mat = Matrix<General, Ref, Ref, AT>;
62 public:
63 Iterator(std::conditional_t<Const, const Mat, Mat> &mat_, int row_, int col_) : mat(mat_), row(row_), col(col_) {}
64 bool operator!=(const Iterator& x) const { return row!=x.row || col!=x.col; }
65 Iterator& operator++() { ++row; if(row>=mat.m) { row=0; ++col; } return *this; }
66 std::conditional_t<Const, const AT, AT>& operator*() { return mat.e(row,col); }
67 std::conditional_t<Const, const AT, AT>& operator->() { return mat.e(row,col); }
68 const AT& operator*() const { return mat.e(row,col); }
69 const AT& operator->() const { return mat.e(row,col); }
70 private:
71 std::conditional_t<Const, const Mat, Mat> &mat;
72 int row;
73 int col;
74 };
75
76 Memory<AT> memory;
77 AT *ele;
78 int m{0};
79 int n{0};
80 int lda{0};
81
82 template <class Type, class Row, class Col> inline Matrix<General,Ref,Ref,AT>& copy(const Matrix<Type,Row,Col,AT> &A);
83 inline Matrix<General,Ref,Ref,AT>& copy(const Matrix<Sparse,Ref,Ref,AT> &A);
84 inline Matrix<General,Ref,Ref,AT>& copy(const Matrix<SymmetricSparse,Ref,Ref,AT> &A);
85
86 const AT* elePtr(int i, int j) const {
87 return ele+i+j*lda;
88 }
89
90 AT* elePtr(int i, int j) {
91 return ele+i+j*lda;
92 }
93
95
96 public:
97
102 explicit Matrix() : memory(), ele(nullptr) { }
103
104 // Works with -std=gnu++0x only
105 // template<class Ini=All<AT>>
106 // Matrix(int m_, int n_, Ini ini=All<AT>()) : memory(m_*n_), ele((AT*)memory.get()), m(m_), n(n_), lda(m_) {
107 // init(ini);
108 // }
109
110 explicit Matrix(int m_, int n_, Noinit) : memory(m_*n_), ele((AT*)memory.get()), m(m_), n(n_), lda(m_) { }
111 explicit Matrix(int m_, int n_, Init ini=INIT, const AT &a=AT()) : memory(m_*n_), ele((AT*)memory.get()), m(m_), n(n_), lda(m_) { init0(a); }
112 explicit Matrix(int m_, int n_, Eye ini, const AT &a=1) : memory(m_*n_), ele((AT*)memory.get()), m(m_), n(n_), lda(m_) { init(ini,a); }
113
119 Matrix(const Matrix<General,Ref,Ref,AT> &A) : memory(A.m*A.n), ele((AT*)memory.get()) ,m(A.m), n(A.n), lda(m) {
120 copy(A);
121 }
122
128 template<class Row, class Col>
129 Matrix(const Matrix<General,Row,Col,AT> &A) : memory(A.rows()*A.cols()), ele((AT*)memory.get()), m(A.rows()), n(A.cols()), lda(m) {
130 copy(A);
131 }
132
138 template<class Type, class Row, class Col>
139 explicit Matrix(const Matrix<Type,Row,Col,AT> &A) : memory(A.rows()*A.cols()), ele((AT*)memory.get()), m(A.rows()), n(A.cols()), lda(m) {
140 copy(A);
141 }
142
150 explicit Matrix(int m_, int n_, AT* ele_) : memory(), ele(ele_), m(m_), n(n_), lda(m_) {
151 }
152
165 Matrix(const std::string &strs);
166 Matrix(const char *strs);
167
170 ~Matrix() = default;
171
172 using iterator = Iterator<false>;
173 using const_iterator = Iterator<true>;
174 iterator begin() { return iterator(*this,0,0); }
175 iterator end() { return iterator(*this,0,n); }
176 const_iterator begin() const { return const_iterator(*this,0,0); }
177 const_iterator end() const { return const_iterator(*this,0,n); }
178
179 Matrix<General,Ref,Ref,AT>& resize(int m_, int n_, Noinit) {
180 m = m_; n = n_; lda = m;
181 memory.resize(m*n);
182 ele = (AT*)memory.get();
183 return *this;
184 }
185
186 Matrix<General,Ref,Ref,AT>& resize(int m, int n, Init ini=INIT, const AT &a=AT()) { return resize(m,n,Noinit()).init0(a); }
187
188 Matrix<General,Ref,Ref,AT>& resize(int m, int n, Eye ini, const AT &a=1) { return resize(m,n,Noinit()).init(ini,a); }
189
196 Matrix<General,Ref,Ref,AT>& operator=(const Matrix<General,Ref,Ref,AT> &A) {
197 FMATVEC_ASSERT(m == A.rows(), AT);
198 FMATVEC_ASSERT(n == A.cols(), AT);
199 return copy(A);
200 }
201
208 template<class Type, class Row, class Col>
209 Matrix<General,Ref,Ref,AT>& operator=(const Matrix<Type,Row,Col,AT> &A) {
210 FMATVEC_ASSERT(m == A.rows(), AT);
211 FMATVEC_ASSERT(n == A.cols(), AT);
212 return copy(A);
213 }
214
221 inline Matrix<General,Ref,Ref,AT>& operator&=(Matrix<General,Ref,Ref,AT> &A) {
222 m=A.m;
223 n=A.n;
224 memory = A.memory;
225 ele = A.ele;
226 lda = A.lda;
227 return *this;
228 }
229
236 template<class Type, class Row, class Col>
237 inline Matrix<General,Ref,Ref,AT>& operator<<=(const Matrix<Type,Row,Col,AT> &A) {
238 if(m!=A.rows() || n!=A.cols()) resize(A.rows(),A.cols(),NONINIT);
239 return copy(A);
240 }
241
242 template <class AT2>
243 operator Matrix<General,Ref,Ref,AT2>() const {
244 Matrix<General,Ref,Ref,AT2> ret(rows(), cols());
245 for(size_t r=0; r<rows(); ++r)
246 for(size_t c=0; c<cols(); ++c)
247 ret(r,c) = (*this)(r,c);
248 return ret;
249 }
250
260 AT& operator()(int i, int j) {
261 FMATVEC_ASSERT(i>=0, AT);
262 FMATVEC_ASSERT(j>=0, AT);
263 FMATVEC_ASSERT(i<m, AT);
264 FMATVEC_ASSERT(j<n, AT);
265
266 return e(i,j);
267 }
268
273 const AT& operator()(int i, int j) const {
274 FMATVEC_ASSERT(i>=0, AT);
275 FMATVEC_ASSERT(j>=0, AT);
276 FMATVEC_ASSERT(i<m, AT);
277 FMATVEC_ASSERT(j<n, AT);
278
279 return e(i,j);
280 }
281
282 AT& e(int i, int j) {
283 return ele[i+j*lda];
284 }
285
286 const AT& e(int i, int j) const {
287 return ele[i+j*lda];
288 }
289
295 AT* operator()() {return ele;}
296
301 const AT* operator()() const {return ele;}
302
307 int rows() const {return m;}
308
313 int cols() const {return n;}
314
319 int ldim() const {return lda;}
320
327 const CBLAS_TRANSPOSE blasTrans() const {
328 return CblasNoTrans;
329 }
330
338 const CBLAS_ORDER blasOrder() const {
339 return CblasColMajor;
340 }
341
373 inline const Matrix<General,Ref,Ref,AT> operator()(const Range<Var,Var> &I, const Range<Var,Var> &J) const;
374
377 inline const SquareMatrix<Ref,AT> operator()(const Range<Var,Var> &I) const;
378
381 inline const RowVector<Ref,AT> row(int i) const;
382
385 inline const Vector<Ref,AT> col(int j) const;
386
394 inline Matrix<General,Ref,Ref,AT>& init(const AT &val=AT());
395 inline Matrix<General,Ref,Ref,AT>& init(Init, const AT &a=AT()) { return init(a); }
396 inline Matrix<General,Ref,Ref,AT>& init(Eye, const AT &val=1);
397 inline Matrix<General,Ref,Ref,AT>& init(Noinit, const AT &a=AT()) { return *this; }
398 inline Matrix<General,Ref,Ref,AT>& init0(const AT &val=AT());
399 inline Matrix<General,Ref,Ref,AT>& init0(Init, const AT &a=AT()) { return init0(a); }
400
405 explicit inline operator std::vector<std::vector<AT>>() const;
406
412 explicit inline Matrix(const std::vector<std::vector<AT>> &m);
413
414 // /*! \brief Cast to AT.
415 // *
416 // * \return The AT representation of the matrix
417 // * */
418 // explicit operator AT() const {
419 // FMATVEC_ASSERT(m==1, AT);
420 // FMATVEC_ASSERT(n==1, AT);
421 // return e(0);
422 // }
423 //
424 // /*! \brief AT Constructor.
425 // * Constructs and initializes a matrix with a AT object.
426 // * \param x The AT the matrix will be initialized with.
427 // * */
428 // explicit Matrix(const AT &x) : memory(1), ele((AT*)memory.get()), m(1), n(1), lda(1) { ele[0] = x; }
429
430 const Matrix<General,Ref,Ref,AT> T() const;
431
432 template<class Row> inline void set(int j, const Vector<Row,AT> &x);
433
434 template<class Col> inline void set(int i, const RowVector<Col,AT> &x);
435
436 template<class Type, class Row, class Col> inline void set(const Range<Var,Var> &I, const Range<Var,Var> &J, const Matrix<Type,Row,Col,AT> &A);
437
438 template<class Row> inline void add(int j, const Vector<Row,AT> &x);
439
440 template<class Col> inline void add(int i, const RowVector<Col,AT> &x);
441
442 template<class Type, class Row, class Col> inline void add(const Range<Var,Var> &I, const Range<Var,Var> &J, const Matrix<Type,Row,Col,AT> &A);
443
444 template<class Type, class Row, class Col> inline void sub(const Range<Var,Var> &I, const Range<Var,Var> &J, const Matrix<Type,Row,Col,AT> &A);
445
446 inline const Matrix<General,Ref,Ref,AT> operator()(const Indices &I, const Indices &J) const;
447
448 template<class Type, class Row, class Col> inline void set(const Indices &I, const Indices &J, const Matrix<Type,Row,Col,AT> &A);
449
450 template<class Row> inline void set(const Indices &I, int j, const Vector<Row,AT> &x);
451
452 inline void ref(Matrix<General,Ref,Ref,AT> &A, const Range<Var,Var> &I, const Range<Var,Var> &J);
453
454 inline void ref(Matrix<Symmetric,Ref,Ref,AT> &A, const Range<Var,Var> &I, const Range<Var,Var> &J);
455
462 int nonZeroElements() const;
463 };
464
465 template <class AT>
466 Matrix<General,Ref,Ref,AT>::Matrix(const std::string &strs) : memory(), ele(nullptr) {
467 std::istringstream iss(strs);
468 iss>>*this;
469
470 // check end of stream
471 iss>>std::ws;
472 if(!iss.eof())
473 throw std::runtime_error("Input not fully read.");
474 }
475 template <class AT> Matrix<General,Ref,Ref,AT>::Matrix(const char *strs) : Matrix<General,Ref,Ref,AT>::Matrix(std::string(strs)) {}
476
477 template <class AT>
478 inline Matrix<General,Ref,Ref,AT>& Matrix<General,Ref,Ref,AT>::init0(const AT &val) {
479 for(int i=0; i<m*n; i++) ele[i]=val;
480 return *this;
481 }
482
483 template <class AT>
484 inline Matrix<General,Ref,Ref,AT>& Matrix<General,Ref,Ref,AT>::init(const AT &val) {
485 for(int i=0; i<rows(); i++)
486 for(int j=0; j<cols(); j++)
487 e(i,j) = val;
488 return *this;
489 }
490
491 template <class AT>
492 inline Matrix<General,Ref,Ref,AT>& Matrix<General,Ref,Ref,AT>::init(Eye, const AT &val) {
493 for(int i=0; i<m; i++)
494 for(int j=0; j<n; j++)
495 e(i,j) = (i==j) ? val : 0;
496 return *this;
497 }
498
499 template <class AT>
500 inline const Matrix<General,Ref,Ref,AT> Matrix<General,Ref,Ref,AT>::operator()(const Range<Var,Var> &I, const Range<Var,Var> &J) const {
501 FMATVEC_ASSERT(I.end()<m, AT);
502 FMATVEC_ASSERT(J.end()<n, AT);
503 Matrix<General,Ref,Ref,AT> A(I.size(),J.size(),NONINIT);
504
505 for(int i=0; i<A.rows(); i++)
506 for(int j=0; j<A.cols(); j++)
507 A.e(i,j) = e(I.start()+i,J.start()+j);
508
509 return A;
510 }
511
512 template <class AT>
513 inline const SquareMatrix<Ref,AT> Matrix<General,Ref,Ref,AT>::operator()(const Range<Var,Var> &I) const {
514 FMATVEC_ASSERT(I.end()<m, AT);
515 SquareMatrix<Ref,AT> A(I.size(),NONINIT);
516
517 for(int i=0; i<A.rows(); i++)
518 for(int j=0; j<A.cols(); j++)
519 A.e(i,j) = e(I.start()+i,I.start()+j);
520
521 return A;
522 }
523
524 template <class AT>
525 inline const RowVector<Ref,AT> Matrix<General,Ref,Ref,AT>::row(int i) const {
526
527 FMATVEC_ASSERT(i>=0, AT);
528 FMATVEC_ASSERT(i<m, AT);
529
530 RowVector<Ref,AT> x(n,NONINIT);
531
532 for(int j=0; j<n; j++)
533 x.e(j) = e(i,j);
534
535 return x;
536 }
537
538 template <class AT>
539 inline const Vector<Ref,AT> Matrix<General,Ref,Ref,AT>::col(int j) const {
540
541 FMATVEC_ASSERT(j>=0, AT);
542 FMATVEC_ASSERT(j<n, AT);
543
544 Vector<Ref,AT> x(m,NONINIT);
545
546 for(int i=0; i<m; i++)
547 x.e(i) = e(i,j);
548
549 return x;
550 }
551
552 template <class AT>
553 inline Matrix<General,Ref,Ref,AT>::operator std::vector<std::vector<AT>>() const {
554 std::vector<std::vector<AT>> ret(rows(),std::vector<AT>(cols()));
555 for(int r=0; r<rows(); r++) {
556 for(int c=0; c<cols(); c++)
557 ret[r][c]= e(r,c);
558 }
559 return ret;
560 }
561
562 template <class AT>
563 inline Matrix<General,Ref,Ref,AT>::Matrix(const std::vector<std::vector<AT>> &m) : memory(m.size()>0?m.size()*m[0].size():0), ele((AT*)memory.get()), m(static_cast<int>(m.size())), n(static_cast<int>(m.size()>0?m[0].size():0)), lda(static_cast<int>(m.size())) {
564 for(int r=0; r<rows(); r++) {
565 if(static_cast<int>(m[r].size())!=cols())
566 throw std::runtime_error("The rows of the input have different length.");
567 for(int c=0; c<cols(); c++)
568 e(r,c)=m[r][c];
569 }
570 }
571
573
574 template <class AT> template <class Type, class Row, class Col>
575 inline Matrix<General,Ref,Ref,AT>& Matrix<General,Ref,Ref,AT>::copy(const Matrix<Type,Row,Col,AT> &A) {
576 for(int i=0; i<m; i++)
577 for(int j=0; j<n; j++)
578 e(i,j) = A.e(i,j);
579 return *this;
580 }
581
582 template <class AT>
583 inline Matrix<General,Ref,Ref,AT>& Matrix<General,Ref,Ref,AT>::copy(const Matrix<Sparse,Ref,Ref,AT> &A) {
584 for(int i=0; i<m; i++) {
585 for(int j=0; j<n; j++)
586 e(i,j) = 0;
587 for(int j=A.Ip()[i]; j<A.Ip()[i+1]; j++)
588 e(i,A.Jp()[j]) = A()[j];
589 }
590 return *this;
591 }
592
593 template <class AT>
594 inline Matrix<General,Ref,Ref,AT>& Matrix<General,Ref,Ref,AT>::copy(const Matrix<SymmetricSparse,Ref,Ref,AT> &A) {
595 for(int i=0; i<m; i++) {
596 for(int j=0; j<n; j++)
597 e(i,j) = 0;
598 e(i,A.Jp()[A.Ip()[i]]) = A()[A.Ip()[i]];
599 }
600 for(int i=0; i<m; i++) {
601 for(int j=A.Ip()[i]+1; j<A.Ip()[i+1]; j++) {
602 e(i,A.Jp()[j]) = A()[j];
603 e(A.Jp()[j],i) = A()[j];
604 }
605 }
606 return *this;
607 }
608
609 template <class AT>
610 inline const Matrix<General,Ref,Ref,AT> Matrix<General,Ref,Ref,AT>::T() const {
611 Matrix<General,Ref,Ref,AT> A(n,m,NONINIT);
612 for(int i=0; i<n; i++)
613 for(int j=0; j<m; j++)
614 A.e(i,j) = e(j,i);
615 return A;
616 }
617
618 template <class AT> template <class Row>
619 inline void Matrix<General,Ref,Ref,AT>::set(int j, const Vector<Row,AT> &x) {
620 FMATVEC_ASSERT(j<cols(), AT);
621 FMATVEC_ASSERT(rows()==x.size(), AT);
622 for(int i=0; i<rows(); i++)
623 e(i,j) = x.e(i);
624 }
625
626 template <class AT> template <class Col>
627 inline void Matrix<General,Ref,Ref,AT>::set(int i, const RowVector<Col,AT> &x) {
628 FMATVEC_ASSERT(i<rows(), AT);
629 FMATVEC_ASSERT(cols()==x.size(), AT);
630 for(int j=0; j<cols(); j++)
631 e(i,j) = x.e(j);
632 }
633
634 template <class AT> template<class Type, class Row, class Col>
635 inline void Matrix<General,Ref,Ref,AT>::set(const Range<Var,Var> &I, const Range<Var,Var> &J, const Matrix<Type,Row,Col,AT> &A) {
636
637 FMATVEC_ASSERT(I.end()<rows(), AT);
638 FMATVEC_ASSERT(J.end()<cols(), AT);
639 FMATVEC_ASSERT(I.size()==A.rows(), AT);
640 FMATVEC_ASSERT(J.size()==A.cols(), AT);
641
642 for(int i=I.start(), ii=0; i<=I.end(); i++, ii++)
643 for(int j=J.start(), jj=0; j<=J.end(); j++, jj++)
644 e(i,j) = A.e(ii,jj);
645 }
646
647 template <class AT> template <class Row>
648 inline void Matrix<General,Ref,Ref,AT>::add(int j, const Vector<Row,AT> &x) {
649 FMATVEC_ASSERT(j<cols(), AT);
650 FMATVEC_ASSERT(rows()==x.size(), AT);
651 for(int i=0; i<rows(); i++)
652 e(i,j) += x.e(i);
653 }
654
655 template <class AT> template <class Col>
656 inline void Matrix<General,Ref,Ref,AT>::add(int i, const RowVector<Col,AT> &x) {
657 FMATVEC_ASSERT(i<rows(), AT);
658 FMATVEC_ASSERT(cols()==x.size(), AT);
659 for(int j=0; j<cols(); j++)
660 e(i,j) += x.e(j);
661 }
662
663 template <class AT> template<class Type, class Row, class Col>
664 inline void Matrix<General,Ref,Ref,AT>::add(const Range<Var,Var> &I, const Range<Var,Var> &J, const Matrix<Type,Row,Col,AT> &A) {
665
666 FMATVEC_ASSERT(I.end()<rows(), AT);
667 FMATVEC_ASSERT(J.end()<cols(), AT);
668 FMATVEC_ASSERT(I.size()==A.rows(), AT);
669 FMATVEC_ASSERT(J.size()==A.cols(), AT);
670
671 for(int i=I.start(), ii=0; i<=I.end(); i++, ii++)
672 for(int j=J.start(), jj=0; j<=J.end(); j++, jj++)
673 e(i,j) += A.e(ii,jj);
674 }
675
676 template <class AT> template<class Type, class Row, class Col>
677 inline void Matrix<General,Ref,Ref,AT>::sub(const Range<Var,Var> &I, const Range<Var,Var> &J, const Matrix<Type,Row,Col,AT> &A) {
678
679 FMATVEC_ASSERT(I.end()<rows(), AT);
680 FMATVEC_ASSERT(J.end()<cols(), AT);
681 FMATVEC_ASSERT(I.size()==A.rows(), AT);
682 FMATVEC_ASSERT(J.size()==A.cols(), AT);
683
684 for(int i=I.start(), ii=0; i<=I.end(); i++, ii++)
685 for(int j=J.start(), jj=0; j<=J.end(); j++, jj++)
686 e(i,j) -= A.e(ii,jj);
687 }
688
689 template <class AT>
690 inline const Matrix<General,Ref,Ref,AT> Matrix<General,Ref,Ref,AT>::operator()(const Indices &I, const Indices &J) const {
691 FMATVEC_ASSERT(I.max()<rows(), AT);
692 FMATVEC_ASSERT(J.max()<cols(), AT);
693
694 Matrix<General,Ref,Ref,AT> A(I.size(),J.size(),NONINIT);
695
696 for(int i=0; i<A.rows(); i++)
697 for(int j=0; j<A.cols(); j++)
698 A.e(i,j) = e(I[i],J[j]);
699
700 return A;
701 }
702
703 template <class AT> template <class Type, class Row, class Col>
704 inline void Matrix<General,Ref,Ref,AT>::set(const Indices &I, const Indices &J, const Matrix<Type,Row,Col,AT> &A) {
705 FMATVEC_ASSERT(I.max()<rows(), AT);
706 FMATVEC_ASSERT(J.max()<cols(), AT);
707 FMATVEC_ASSERT(I.size()==A.rows(), AT);
708 FMATVEC_ASSERT(J.size()==A.cols(), AT);
709 for(int i=0; i<I.size(); i++)
710 for(int j=0; j<J.size(); j++)
711 e(I[i],J[j]) = A.e(i,j);
712 }
713
714 template <class AT> template<class Row>
715 inline void Matrix<General,Ref,Ref,AT>::set(const Indices &I, int j, const Vector<Row,AT> &x) {
716 FMATVEC_ASSERT(I.max()<rows(), AT);
717 FMATVEC_ASSERT(j<cols(), AT);
718 FMATVEC_ASSERT(I.size()==x.size(), AT);
719 for(int i=0; i<I.size(); i++)
720 e(I[i],j) = x.e(i);
721 }
722
723 template <class AT>
724 inline void Matrix<General,Ref,Ref,AT>::ref(Matrix<General,Ref,Ref,AT> &A, const Range<Var,Var> &I, const Range<Var,Var> &J) {
725 FMATVEC_ASSERT(I.end()<A.rows(), AT);
726 FMATVEC_ASSERT(J.end()<A.cols(), AT);
727 m=I.size();
728 n=J.size();
729 memory = A.memory;
730 ele = A.elePtr(I.start(),J.start());
731 lda = A.lda;
732 }
733
734 template <class AT>
735 inline void Matrix<General,Ref,Ref,AT>::ref(Matrix<Symmetric,Ref,Ref,AT> &A, const Range<Var,Var> &I, const Range<Var,Var> &J) {
736 FMATVEC_ASSERT(I.end()<A.size(), AT);
737 FMATVEC_ASSERT(J.end()<A.size(), AT);
738 FMATVEC_ASSERT(I.start()>=J.start(), AT);
739 FMATVEC_ASSERT(J.end()<=I.start(), AT);
740 m=I.size();
741 n=J.size();
742 memory = A.memory;
743 ele = A.elePtr(I.start(),J.start());
744 lda = A.lda;
745 }
746
747 template <class AT>
748 int Matrix<General,Ref,Ref,AT>::nonZeroElements() const {
749 int k = 0;
750 for (int i = 0; i < m; i++) {
751 for (int j = 0; j < n; j++) {
752 if (fabs(e(i, j)) > 1e-16 || i==j) {
753 k++;
754 }
755 }
756 }
757 return k;
758 }
759
761
762 }
763
764 #endif
fmatvec::Range< Var, Var >::end
int end() const
Last element.
Definition: range.h:91
fmatvec::Matrix< General, Ref, Ref, AT >::blasOrder
const CBLAS_ORDER blasOrder() const
Storage convention.
Definition: general_matrix.h:338
fmatvec::Matrix< General, Ref, Ref, AT >::Matrix
Matrix(const Matrix< General, Row, Col, AT > &A)
Copy Constructor.
Definition: general_matrix.h:129
fmatvec::Matrix
This is the basic matrix class for arbitrary matrices.
Definition: fmatvec.h:42
fmatvec::Symmetric
Shape class for symmetric matrices.
Definition: types.h:128
fmatvec::Matrix< General, Ref, Ref, AT >::rows
int rows() const
Number of rows.
Definition: general_matrix.h:307
fmatvec::Matrix< SymmetricSparse, Ref, Ref, AT >
This is a matrix class for sparse quadratic matrices.
Definition: symmetric_sparse_matrix.h:44
fmatvec::Matrix::operator()
AT & operator()(int i, int j)
Standard constructor.
Definition: matrix.h:84
fmatvec::Init
Definition: types.h:89
fmatvec::Matrix< SymmetricSparse, Ref, Ref, AT >::rows
int rows() const
Number of rows.
Definition: symmetric_sparse_matrix.h:270
fmatvec::Matrix< General, Ref, Ref, AT >::operator()
AT & operator()(int i, int j)
Element operator.
Definition: general_matrix.h:260
fmatvec::Matrix< SymmetricSparse, Ref, Ref, AT >::cols
int cols() const
Number of columns.
Definition: symmetric_sparse_matrix.h:276
fmatvec::Matrix< General, Ref, Ref, AT >::Matrix
Matrix(const Matrix< General, Ref, Ref, AT > &A)
Copy Constructor.
Definition: general_matrix.h:119
fmatvec::RowVector
Definition: fmatvec.h:51
fmatvec::Matrix< General, Ref, Ref, AT >::blasTrans
const CBLAS_TRANSPOSE blasTrans() const
Transposed status.
Definition: general_matrix.h:327
fmatvec::Matrix< Sparse, Ref, Ref, AT >
This is a matrix class for sparse quadratic matrices.
Definition: sparse_matrix.h:44
fmatvec::Ref
Definition: types.h:98
fmatvec::Matrix< General, Ref, Ref, AT >::T
const Matrix< General, Ref, Ref, AT > T() const
Cast to AT.
fmatvec::RowVector< Ref, AT >
This is a rowvector class of general shape in dense storage format.
Definition: row_vector.h:36
fmatvec::Range< Var, Var >::start
int start() const
First element.
Definition: range.h:85
fmatvec::SquareMatrix
Definition: fmatvec.h:45
fmatvec::Vector
Definition: fmatvec.h:48
fmatvec::Matrix< General, Ref, Ref, AT >::operator()
const AT * operator()() const
Pointer operator.
Definition: general_matrix.h:301
fmatvec::Matrix::cols
int cols() const
Number of columns.
fmatvec::Matrix< General, Ref, Ref, AT >
This is a matrix class for general matrices.
Definition: general_matrix.h:42
fmatvec
Namespace fmatvec.
Definition: _memory.cc:28
fmatvec::Matrix< General, Ref, Ref, AT >::cols
int cols() const
Number of columns.
Definition: general_matrix.h:313
fmatvec::Matrix< General, Ref, Ref, AT >::nonZeroElements
int nonZeroElements() const
Count nonzero elements.
fmatvec::operator!=
bool operator!=(const Matrix< Type1, Row1, Col1, AT > &A, const Matrix< Type2, Row2, Col2, AT > &B)
Matrix-matrix comparison.
Definition: linear_algebra.h:1935
fmatvec::Range< Var, Var >::size
int size() const
Size.
Definition: range.h:97
fmatvec::Matrix< General, Ref, Ref, AT >::operator=
Matrix< General, Ref, Ref, AT > & operator=(const Matrix< Type, Row, Col, AT > &A)
Assignment operator.
Definition: general_matrix.h:209
fmatvec::Matrix< General, Ref, Ref, AT >::operator()
AT * operator()()
Pointer operator.
Definition: general_matrix.h:295
fmatvec::Matrix< General, Ref, Ref, AT >::operator&=
Matrix< General, Ref, Ref, AT > & operator&=(Matrix< General, Ref, Ref, AT > &A)
Reference operator.
Definition: general_matrix.h:221
fmatvec::Matrix< General, Ref, Ref, AT >::ldim
int ldim() const
Leading dimension.
Definition: general_matrix.h:319
fmatvec::Matrix< General, Ref, Ref, AT >::Matrix
Matrix(int m_, int n_, AT *ele_)
Regular Constructor.
Definition: general_matrix.h:150
fmatvec::Noinit
Definition: types.h:88
fmatvec::General
Shape class for general matrices.
Definition: types.h:112
fmatvec::Matrix< General, Ref, Ref, AT >::operator()
const AT & operator()(int i, int j) const
Element operator.
Definition: general_matrix.h:273
fmatvec::Vector< Ref, AT >
This is a vector class of general shape in dense storage format.
Definition: vector.h:39
fmatvec::SquareMatrix< Ref, AT >
This is a matrix class of general quadratic matrices.
Definition: square_matrix.h:39
fmatvec::Matrix< General, Ref, Ref, AT >::Matrix
Matrix()
Standard constructor.
Definition: general_matrix.h:102
fmatvec::Matrix< General, Ref, Ref, AT >::operator=
Matrix< General, Ref, Ref, AT > & operator=(const Matrix< General, Ref, Ref, AT > &A)
Assignment operator.
Definition: general_matrix.h:196
fmatvec::Matrix::rows
int rows() const
Number of rows.
fmatvec::Matrix< General, Ref, Ref, AT >::Matrix
Matrix(const Matrix< Type, Row, Col, AT > &A)
Copy Constructor.
Definition: general_matrix.h:139
fmatvec::Range< Var, Var >
This is an index class for creating submatrices.
Definition: range.h:44
fmatvec::Eye
Definition: types.h:90