All Classes Namespaces Functions Variables Typedefs Enumerations Pages
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 "matrix.h"
28 #include <stdlib.h>
29 
30 namespace fmatvec {
31 
40  template <class AT> class Matrix<General,Ref,Ref,AT> {
41 
42  public:
43 
44  typedef AT AtomicType;
45 
47 
48  friend class Vector<Ref,AT>;
49  friend class RowVector<Ref,AT>;
50  friend class SquareMatrix<Ref,AT>;
51  friend class Matrix<Symmetric,Ref,Ref,AT>;
52 
53  template <class T> friend Matrix<General,Ref,Ref,T> trans(const Matrix<General,Ref,Ref,T> &A);
54 
55  protected:
56 
57  Memory<AT> memory;
58  AT *ele;
59  int m;
60  int n;
61  int lda;
62  bool tp;
63 
64  template <class Type, class Row, class Col> inline void deepCopy(const Matrix<Type,Row,Col,AT> &A);
65  inline void deepCopy(const Matrix<General,Ref,Ref,AT> &A);
66  inline void deepCopy(const Matrix<Symmetric,Ref,Ref,AT> &A);
67 
68  const AT* elePtr(int i, int j) const {
69  return tp ? ele+i*lda+j : ele+i+j*lda;
70  };
71 
72  AT* elePtr(int i, int j) {
73  return tp ? ele+i*lda+j : ele+i+j*lda;
74  };
75 
76  Matrix(int m_, int n_, int lda_, bool tp_, Memory<AT> memory_, const AT* ele_) : memory(memory_), ele((AT*)ele_), m(m_), n(n_), lda(lda_), tp(tp_) {
77  }
78 
80 
81  public:
82 
87  Matrix() : memory(), ele(0), m(0), n(0), lda(0), tp(false) { }
88 
89 // Works with -std=gnu++0x only
90 // template<class Ini=All<AT> >
91 // Matrix(int m_, int n_, Ini ini=All<AT>()) : memory(m_*n_), ele((AT*)memory.get()), m(m_), n(n_), lda(m_), tp(false) {
92 // init(ini);
93 // }
94 
95  Matrix(int m_, int n_, Noinit) : memory(m_*n_), ele((AT*)memory.get()), m(m_), n(n_), lda(m_), tp(false) { }
96  Matrix(int m_, int n_, Init ini=INIT, const AT &a=0) : memory(m_*n_), ele((AT*)memory.get()), m(m_), n(n_), lda(m_), tp(false) { init0(a); }
97  Matrix(int m_, int n_, Eye ini, const AT &a=1) : memory(m_*n_), ele((AT*)memory.get()), m(m_), n(n_), lda(m_), tp(false) { init(ini,a); }
98 
106  Matrix(const Matrix<General,Ref,Ref,AT> &A) : memory(A.memory), ele(A.ele) ,m(A.m), n(A.n), lda(A.lda), tp(A.tp) {
107  }
108 
116  Matrix(int m_, int n_, AT* ele_) : memory(), ele(ele_), m(m_), n(n_), lda(m_), tp(false) {
117  }
118 
131  Matrix(const char *str);
132 
136  }
137 
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), tp(false) {
140 
141  deepCopy(A);
142  }
143 
144  Matrix<General,Ref,Ref,AT>& resize() {
145  m = n = lda = 0;
146  tp = false;
147  memory.resize(0);
148  ele = 0;
149  return *this;
150  }
151 
152  Matrix<General,Ref,Ref,AT>& resize(int m_, int n_, Noinit) {
153  m = m_; n = n_; lda = m;
154  tp = false;
155  memory.resize(m*n);
156  ele = (AT*)memory.get();
157  return *this;
158  }
159 
160  Matrix<General,Ref,Ref,AT>& resize(int m, int n, Init ini=INIT, const AT &a=0) { return resize(m,n,Noinit()).init0(a); }
161 
162  Matrix<General,Ref,Ref,AT>& resize(int m, int n, Eye ini, const AT &a=1) { return resize(m,n,Noinit()).init(ini,a); }
163 
172  inline Matrix<General,Ref,Ref,AT>& operator=(const Matrix<General,Ref,Ref,AT> &A);
173 
174  template<class Type, class Row, class Col>
175  Matrix<General,Ref,Ref,AT>& operator=(const Matrix<Type,Row,Col,AT> &A);
176 
183  template<class Type, class Row, class Col>
184  inline Matrix<General,Ref,Ref,AT>& operator<<(const Matrix<Type,Row,Col,AT> &A);
185 
192  inline Matrix<General,Ref,Ref,AT>& operator>>(const Matrix<General,Ref,Ref,AT> &A);
193 
205  AT& operator()(int i, int j) {
206 #ifndef FMATVEC_NO_BOUNDS_CHECK
207  assert(i>=0);
208  assert(j>=0);
209  assert(i<m);
210  assert(j<n);
211 #endif
212 
213  return e(i,j);
214  };
215 
220  const AT& operator()(int i, int j) const {
221 #ifndef FMATVEC_NO_BOUNDS_CHECK
222  assert(i>=0);
223  assert(j>=0);
224  assert(i<m);
225  assert(j<n);
226 #endif
227 
228  return e(i,j);
229  };
230 
231  AT& er(int i, int j) {
232  return ele[i+j*lda];
233  };
234 
235  const AT& er(int i, int j) const {
236  return ele[i+j*lda];
237  };
238 
239  AT& et(int i, int j) {
240  return ele[i*lda+j];
241  };
242 
243  const AT& et(int i, int j) const {
244  return ele[i*lda+j];
245  };
246 
247  AT& e(int i, int j) {
248  return tp ? et(i,j) : er(i,j);
249  };
250 
251  const AT& e(int i, int j) const {
252  return tp ? et(i,j) : er(i,j);
253  };
254 
260  AT* operator()() {return ele;};
261 
266  const AT* operator()() const {return ele;};
267 
272  int rows() const {return m;};
273 
278  int cols() const {return n;};
279 
284  int ldim() const {return lda;};
285 
291  bool transposed() const {return tp;};
292 
299  const CBLAS_TRANSPOSE blasTrans() const {
300  return (tp)? CblasTrans : CblasNoTrans;
301  };
302 
310  const CBLAS_ORDER blasOrder() const {
311  return CblasColMajor;
312  };
313 
325  inline Matrix<General,Ref,Ref,AT> operator()(int i1, int j1, int i2, int j2);
326 
331  inline const Matrix<General,Ref,Ref,AT> operator()(int i1, int j1, int i2, int j2) const;
332 
360 
365  inline const Matrix<General,Ref,Ref,AT> operator()(const Range<Var,Var> &I, const Range<Var,Var> &J) const;
366 
376 
381  inline const SquareMatrix<Ref,AT> operator()(const Range<Var,Var> &I) const;
382 
390  inline Vector<Ref,AT> col(int i);
391 
396  inline const Vector<Ref,AT> col(int i) const;
397 
405  inline RowVector<Ref,AT> row(int i);
406 
411  inline const RowVector<Ref,AT> row(int i) const;
412 
418  inline Matrix<General,Ref,Ref,AT> copy() const;
419 
427  inline Matrix<General,Ref,Ref,AT>& init(const AT &a=0);
428  inline Matrix<General,Ref,Ref,AT>& init(Init, const AT &a=0) { return init(a); }
429  inline Matrix<General,Ref,Ref,AT>& init(Eye, const AT &a=1);
430  inline Matrix<General,Ref,Ref,AT>& init(Noinit, const AT &a=0) { return *this; }
431  inline Matrix<General,Ref,Ref,AT>& init0(const AT &a=0);
432  inline Matrix<General,Ref,Ref,AT>& init0(Init, const AT &a=0) { return init0(a); }
433 
438  inline operator std::vector<std::vector<AT> >();
439 
445  inline Matrix(std::vector<std::vector<AT> > m);
446 
447  Matrix<General,Ref,Ref,AT> T() {
448  return Matrix<General,Ref,Ref,AT>(n,m,lda,tp?false:true,memory,ele);
449  };
450 
451  const Matrix<General,Ref,Ref,AT> T() const {
452  return Matrix<General,Ref,Ref,AT>(n,m,lda,tp?false:true,memory,ele);
453  }
454  };
455 
456  template <class AT>
457  Matrix<General,Ref,Ref,AT>::Matrix(const char *strs) : memory(), ele(0), m(0), n(0), lda(0) {
458  // if 'strs' is a single scalar, surround it first with '[' and ']'.
459  // This is more Matlab-like, because e.g. '5' and '[5]' is just the same.
460  // (This functionallitiy is needed e.g. by MBXMLUtils (OpenMBV,MBSim))
461  std::istringstream iss(strs);
462  char c;
463  iss>>c;
464  if(c=='[') iss.str(strs);
465  else iss.str(std::string("[")+strs+"]");
466 
467  int buf=0;
468  iss >> c;
469  iss >> c;
470  if(c!=']') {
471  iss.putback(c);
472  AT x;
473  do {
474  iss >> x;
475  iss >> c;
476  if(c==';') {
477  if(buf)
478  assert(buf == n);
479 
480  buf=n;
481  n=0;
482  m++;
483  }
484  else if(c==',')
485  n++;
486  c='0';
487  } while(iss);
488 
489  n++; m++;
490  lda=m;
491  tp = false;
492  memory.resize(m*n);
493  ele = (AT*)memory.get();
494  iss.clear();
495  iss.seekg(0);
496  iss >> c;
497  for(int i=0; i<m; i++)
498  for(int j=0; j<n; j++) {
499  iss >> er(i,j);
500  iss >> c;
501  }
502  }
503  }
504 
505  template <class AT>
507 
508  m=A.m;
509  n=A.n;
510  memory = A.memory;
511  ele = A.ele;
512  lda = A.lda;
513  tp = A.tp;
514 
515  return *this;
516  }
517 
518  template <class AT>
520 
521  if(!ele) {
522  m = A.rows();
523  n = A.cols();
524  lda = m;
525  tp = false;
526  memory.resize(m*n);
527  ele = (AT*)memory.get();
528  } else {
529 #ifndef FMATVEC_NO_SIZE_CHECK
530  assert(m == A.rows());
531  assert(n == A.cols());
532 #endif
533  }
534 
535  deepCopy(A);
536 
537  return *this;
538  }
539 
540  template <class AT> template< class Type, class Row, class Col>
542 
543  if(!ele) {
544  m = A.rows();
545  n = A.cols();
546  lda = m;
547  tp = false;
548  memory.resize(m*n);
549  ele = (AT*)memory.get();
550  } else {
551 #ifndef FMATVEC_NO_SIZE_CHECK
552  assert(m == A.rows());
553  assert(n == A.cols());
554 #endif
555  }
556 
557  deepCopy(A);
558 
559  return *this;
560  }
561 
562  template <class AT> template< class Type, class Row, class Col>
564 
565  if(m!=A.rows() || n!=A.cols()) {
566  m = A.rows();
567  n = A.cols();
568  lda = m;
569  tp = false;
570  memory.resize(m*n);
571  ele = (AT*)memory.get();
572  }
573 
574  deepCopy(A);
575 
576  return *this;
577  }
578 
579  template <class AT>
581  for(int i=0; i<m*n; i++) ele[i]=val;
582  return *this;
583  }
584 
585  template <class AT>
587 
588  if(tp) {
589  for(int i=0; i<rows(); i++)
590  for(int j=0; j<cols(); j++)
591  et(i,j) = val;
592  }
593  else {
594  for(int i=0; i<rows(); i++)
595  for(int j=0; j<cols(); j++)
596  er(i,j) = val;
597  }
598 
599  return *this;
600  }
601 
602  template <class AT>
604  if(tp)
605  for(int i=0; i<m; i++)
606  for(int j=0; j<n; j++)
607  et(i,j) = (i==j) ? val : 0;
608  else
609  for(int i=0; i<m; i++)
610  for(int j=0; j<n; j++)
611  er(i,j) = (i==j) ? val : 0;
612  return *this;
613  }
614 
615  template <class AT>
617  return operator()(Range<Var,Var>(i1,i2),Range<Var,Var>(j1,j2));
618  }
619 
620  template <class AT>
621  inline const Matrix<General,Ref,Ref,AT> Matrix<General,Ref,Ref,AT>::operator()(int i1, int j1, int i2, int j2) const {
622  return operator()(Range<Var,Var>(i1,i2),Range<Var,Var>(j1,j2));
623  }
624 
625  template <class AT>
627 #ifndef FMATVEC_NO_BOUNDS_CHECK
628  assert(I.end()<m);
629  assert(J.end()<n);
630 #endif
631  return Matrix<General,Ref,Ref,AT>(I.end()-I.start()+1,J.end()-J.start()+1,lda,tp,memory,elePtr(I.start(),J.start()));
632  }
633 
634  template <class AT>
636 #ifndef FMATVEC_NO_BOUNDS_CHECK
637  assert(I.end()<m);
638  assert(J.end()<n);
639 #endif
640  return Matrix<General,Ref,Ref,AT>(I.end()-I.start()+1,J.end()-J.start()+1,lda,tp,memory,elePtr(I.start(),J.start()));
641  }
642 
643  template <class AT>
645 #ifndef FMATVEC_NO_BOUNDS_CHECK
646  assert(I.end()<m);
647 #endif
648  return SquareMatrix<Ref,AT>(I.end()-I.start()+1,lda,tp,memory,elePtr(I.start(),I.start()));
649  }
650 
651  template <class AT>
653 #ifndef FMATVEC_NO_BOUNDS_CHECK
654  assert(I.end()<m);
655 #endif
656  return SquareMatrix<Ref,AT>(I.end()-I.start()+1,lda,tp,memory,elePtr(I.start(),I.start()));
657  }
658 
659  template <class AT>
661 
662 #ifndef FMATVEC_NO_BOUNDS_CHECK
663  assert(i>=0);
664  assert(i<n);
665 #endif
666 
667  return Vector<Ref,AT>(m,lda,tp,memory,elePtr(0,i));
668  }
669 
670  template <class AT>
672 
673 #ifndef FMATVEC_NO_BOUNDS_CHECK
674  assert(i>=0);
675  assert(i<n);
676 #endif
677 
678  return Vector<Ref,AT>(m,lda,tp,memory,elePtr(0,i));
679  }
680 
681  template <class AT>
683 
684 #ifndef FMATVEC_NO_BOUNDS_CHECK
685  assert(i>=0);
686  assert(i<m);
687 #endif
688 
689  return RowVector<Ref,AT>(n,lda,tp,memory,elePtr(i,0));
690  }
691 
692  template <class AT>
694 
695 #ifndef FMATVEC_NO_BOUNDS_CHECK
696  assert(i>=0);
697  assert(i<m);
698 #endif
699 
700  return RowVector<Ref,AT>(n,lda,tp,memory,elePtr(i,0));
701  }
702 
703  template <class AT>
705 
706  Matrix<General,Ref,Ref,AT> A(m,n,NONINIT);
707  A.deepCopy(*this);
708 
709  return A;
710  }
711 
712  template <class AT>
713  inline Matrix<General,Ref,Ref,AT>::operator std::vector<std::vector<AT> >() {
714  std::vector<std::vector<AT> > ret(rows());
715  if(tp) {
716  for(int r=0; r<rows(); r++) {
717  ret[r].resize(cols());
718  for(int c=0; c<cols(); c++)
719  ret[r][c]= et(r,c);
720  }
721  }
722  else {
723  for(int r=0; r<rows(); r++) {
724  ret[r].resize(cols());
725  for(int c=0; c<cols(); c++)
726  ret[r][c]= er(r,c);
727  }
728  }
729  return ret;
730  }
731 
732  template <class AT>
733  inline Matrix<General,Ref,Ref,AT>::Matrix(std::vector<std::vector<AT> > m) : memory(m.size()*m[0].size()), ele((AT*)memory.get()), m(m.size()), n(m[0].size()), lda(m.size()), tp(false) {
734 #ifndef FMATVEC_NO_INITIALIZATION
735  init(0);
736 #endif
737  for(int r=0; r<rows(); r++) {
738  assert(m[r].size()==cols());
739  for(int c=0; c<cols(); c++)
740  er(r,c)=m[r][c];
741  }
742  }
743 
745 
746  template <class AT> template <class Type, class Row, class Col>
748  if(tp) {
749  for(int i=0; i<m; i++)
750  for(int j=0; j<n; j++)
751  et(i,j) = A.e(i,j);
752  }
753  else {
754  for(int i=0; i<m; i++)
755  for(int j=0; j<n; j++)
756  er(i,j) = A.e(i,j);
757  }
758  }
759 
760  template <class AT>
761  inline void Matrix<General,Ref,Ref,AT>::deepCopy(const Matrix<General,Ref,Ref,AT> &A) {
762  if(A.tp) {
763  if(tp) {
764  for(int i=0; i<m; i++)
765  for(int j=0; j<n; j++)
766  et(i,j) = A.et(i,j);
767  }
768  else {
769  for(int i=0; i<m; i++)
770  for(int j=0; j<n; j++)
771  er(i,j) = A.et(i,j);
772  }
773  } else {
774  if(tp) {
775  for(int i=0; i<m; i++)
776  for(int j=0; j<n; j++)
777  et(i,j) = A.er(i,j);
778  }
779  else {
780  for(int i=0; i<m; i++)
781  for(int j=0; j<n; j++)
782  er(i,j) = A.er(i,j);
783  }
784  }
785  }
786 
787  template<class AT>
788  inline void Matrix<General,Ref,Ref,AT>::deepCopy(const Matrix<Symmetric,Ref,Ref,AT> &A) {
789  for(int i=0; i<A.size(); i++) {
790  er(i,i) = A.ej(i,i);
791  for(int j=i+1; j<A.size(); j++)
792  er(i,j) = et(i,j) = A.ej(i,j);
793  }
794  }
795 
797 
798 }
799 
800 #endif
This is the basic matrix class for arbitrary matrices.
Definition: fmatvec.h:41
int ldim() const
Leading dimension.
Definition: general_matrix.h:284
AT & operator()(int i, int j)
Standard constructor.
Definition: matrix.h:86
int rows() const
Number of rows.
Definition: general_matrix.h:272
Definition: matrix.h:39
AT & operator()(int i, int j)
Element operator.
Definition: general_matrix.h:205
Matrix(const Matrix< General, Ref, Ref, AT > &A)
Copy Constructor.
Definition: general_matrix.h:106
Definition: fmatvec.h:50
const CBLAS_TRANSPOSE blasTrans() const
Transposed status.
Definition: general_matrix.h:299
int rows() const
Number of rows.
int cols() const
Number of columns.
Definition: types.h:62
This is a rowvector class of general shape in dense storage format.
Definition: row_vector.h:36
Definition: fmatvec.h:44
Definition: fmatvec.h:47
const CBLAS_ORDER blasOrder() const
Storage convention.
Definition: general_matrix.h:310
This is a matrix class for general matrices.
Definition: general_matrix.h:40
const AT & operator()(int i, int j) const
Element operator.
Definition: general_matrix.h:220
std::istream & operator>>(std::istream &is, Matrix< Type, Row, Col, AT > &A)
Matrix input.
Definition: matrix.h:171
AT * operator()()
Pointer operator.
Definition: general_matrix.h:260
Matrix(int m_, int n_, AT *ele_)
Regular Constructor.
Definition: general_matrix.h:116
Definition: matrix.h:38
RowVector< Ref, AT > trans(const Vector< Ref, AT > &x)
Transpose of a vector.
Definition: linear_algebra.h:1470
This is a matrix class for symmetric matrices.
Definition: symmetric_matrix.h:39
bool transposed() const
Transposed status.
Definition: general_matrix.h:291
Basic shape class for matrices.
Definition: types.h:100
This is a vector class of general shape in dense storage format.
Definition: vector.h:39
This is a matrix class of general quadratic matrices.
Definition: square_matrix.h:39
Matrix()
Standard constructor.
Definition: general_matrix.h:87
int end() const
Last element.
Definition: range.h:95
~Matrix()
Destructor.
Definition: general_matrix.h:135
This is an index class for creating submatrices.
Definition: range.h:44
const AT * operator()() const
Pointer operator.
Definition: general_matrix.h:266
int start() const
First element.
Definition: range.h:89
int cols() const
Number of columns.
Definition: general_matrix.h:278
Definition: matrix.h:40

Impressum / Disclaimer / Datenschutz Generated by doxygen 1.8.5 Valid HTML