All Classes Namespaces Functions 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 
45 
46  friend class Vector<Ref,AT>;
47  friend class RowVector<Ref,AT>;
48  friend class SquareMatrix<Ref,AT>;
49  friend class Matrix<Symmetric,Ref,Ref,AT>;
50 
51  template <class T> friend Matrix<General,Ref,Ref,T> trans(const Matrix<General,Ref,Ref,T> &A);
52 
53  protected:
54 
55  Memory<AT> memory;
56  AT *ele;
57  int m;
58  int n;
59  int lda;
60  bool tp;
61 
62  template <class Type, class Row, class Col> inline void deepCopy(const Matrix<Type,Row,Col,AT> &A);
63  inline void deepCopy(const Matrix<General,Ref,Ref,AT> &A);
64  inline void deepCopy(const Matrix<Symmetric,Ref,Ref,AT> &A);
65 
66  const AT* elePtr(int i, int j) const {
67  return tp ? ele+i*lda+j : ele+i+j*lda;
68  };
69 
70  AT* elePtr(int i, int j) {
71  return tp ? ele+i*lda+j : ele+i+j*lda;
72  };
73 
74  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_) {
75  }
76 
78 
79  public:
80 
85  Matrix() : memory(), ele(0), m(0), n(0), lda(0), tp(false) { }
86 
87 // Works with -std=gnu++0x only
88 // template<class Ini=All<AT> >
89 // Matrix(int m_, int n_, Ini ini=All<AT>()) : memory(m_*n_), ele((AT*)memory.get()), m(m_), n(n_), lda(m_), tp(false) {
90 // init(ini);
91 // }
92 
93  Matrix(int m_, int n_, Noinit) : memory(m_*n_), ele((AT*)memory.get()), m(m_), n(n_), lda(m_), tp(false) { }
94  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); }
95  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); }
96 
104  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) {
105  }
106 
114  Matrix(int m_, int n_, AT* ele_) : memory(), ele(ele_), m(m_), n(n_), lda(m_), tp(false) {
115  }
116 
129  Matrix(const char *str);
130 
134  }
135 
136  template<class Type, class Row, class Col>
137  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) {
138 
139  deepCopy(A);
140  }
141 
142  Matrix<General,Ref,Ref,AT>& resize() {
143  m = n = lda = 0;
144  tp = false;
145  memory.resize(0);
146  ele = 0;
147  return *this;
148  }
149 
150  Matrix<General,Ref,Ref,AT>& resize(int m_, int n_, Noinit) {
151  m = m_; n = n_; lda = m;
152  tp = false;
153  memory.resize(m*n);
154  ele = (AT*)memory.get();
155  return *this;
156  }
157 
158  Matrix<General,Ref,Ref,AT>& resize(int m, int n, Init ini=INIT, const AT &a=0) { return resize(m,n,Noinit()).init0(a); }
159 
160  Matrix<General,Ref,Ref,AT>& resize(int m, int n, Eye ini, const AT &a=1) { return resize(m,n,Noinit()).init(ini,a); }
161 
170  inline Matrix<General,Ref,Ref,AT>& operator=(const Matrix<General,Ref,Ref,AT> &A);
171 
172  template<class Type, class Row, class Col>
173  Matrix<General,Ref,Ref,AT>& operator=(const Matrix<Type,Row,Col,AT> &A);
174 
181  template<class Type, class Row, class Col>
182  inline Matrix<General,Ref,Ref,AT>& operator<<(const Matrix<Type,Row,Col,AT> &A);
183 
190  inline Matrix<General,Ref,Ref,AT>& operator>>(const Matrix<General,Ref,Ref,AT> &A);
191 
203  AT& operator()(int i, int j) {
204 #ifndef FMATVEC_NO_BOUNDS_CHECK
205  assert(i>=0);
206  assert(j>=0);
207  assert(i<m);
208  assert(j<n);
209 #endif
210 
211  return e(i,j);
212  };
213 
218  const AT& operator()(int i, int j) const {
219 #ifndef FMATVEC_NO_BOUNDS_CHECK
220  assert(i>=0);
221  assert(j>=0);
222  assert(i<m);
223  assert(j<n);
224 #endif
225 
226  return e(i,j);
227  };
228 
229  AT& er(int i, int j) {
230  return ele[i+j*lda];
231  };
232 
233  const AT& er(int i, int j) const {
234  return ele[i+j*lda];
235  };
236 
237  AT& et(int i, int j) {
238  return ele[i*lda+j];
239  };
240 
241  const AT& et(int i, int j) const {
242  return ele[i*lda+j];
243  };
244 
245  AT& e(int i, int j) {
246  return tp ? et(i,j) : er(i,j);
247  };
248 
249  const AT& e(int i, int j) const {
250  return tp ? et(i,j) : er(i,j);
251  };
252 
258  AT* operator()() {return ele;};
259 
264  const AT* operator()() const {return ele;};
265 
270  int rows() const {return m;};
271 
276  int cols() const {return n;};
277 
282  int ldim() const {return lda;};
283 
289  bool transposed() const {return tp;};
290 
297  const CBLAS_TRANSPOSE blasTrans() const {
298  return (tp)? CblasTrans : CblasNoTrans;
299  };
300 
308  const CBLAS_ORDER blasOrder() const {
309  return CblasColMajor;
310  };
311 
323  inline Matrix<General,Ref,Ref,AT> operator()(int i1, int j1, int i2, int j2);
324 
329  inline const Matrix<General,Ref,Ref,AT> operator()(int i1, int j1, int i2, int j2) const;
330 
357  inline Matrix<General,Ref,Ref,AT> operator()(const Index &I, const Index &J);
358 
363  inline const Matrix<General,Ref,Ref,AT> operator()(const Index &I, const Index &J) const;
364 
373  inline SquareMatrix<Ref,AT> operator() (const Index &I);
374 
379  inline const SquareMatrix<Ref,AT> operator()(const Index &I) const;
380 
388  inline Vector<Ref,AT> col(int i);
389 
394  inline const Vector<Ref,AT> col(int i) const;
395 
403  inline RowVector<Ref,AT> row(int i);
404 
409  inline const RowVector<Ref,AT> row(int i) const;
410 
416  inline Matrix<General,Ref,Ref,AT> copy() const;
417 
425  inline Matrix<General,Ref,Ref,AT>& init(const AT &a=0);
426  inline Matrix<General,Ref,Ref,AT>& init(Init, const AT &a=0) { return init(a); }
427  inline Matrix<General,Ref,Ref,AT>& init(Eye, const AT &a=1);
428  inline Matrix<General,Ref,Ref,AT>& init(Noinit, const AT &a=0) { return *this; }
429  inline Matrix<General,Ref,Ref,AT>& init0(const AT &a=0);
430  inline Matrix<General,Ref,Ref,AT>& init0(Init, const AT &a=0) { return init0(a); }
431 
436  inline operator std::vector<std::vector<AT> >();
437 
443  inline Matrix(std::vector<std::vector<AT> > m);
444 
445  Matrix<General,Ref,Ref,AT> T() {
446  return Matrix<General,Ref,Ref,AT>(n,m,lda,tp?false:true,memory,ele);
447  };
448 
449  const Matrix<General,Ref,Ref,AT> T() const {
450  return Matrix<General,Ref,Ref,AT>(n,m,lda,tp?false:true,memory,ele);
451  }
452  };
453 
454  template <class AT>
455  Matrix<General,Ref,Ref,AT>::Matrix(const char *strs) : memory(), ele(0), m(0), n(0), lda(0) {
456  // if 'strs' is a single scalar, surround it first with '[' and ']'.
457  // This is more Matlab-like, because e.g. '5' and '[5]' is just the same.
458  // (This functionallitiy is needed e.g. by MBXMLUtils (OpenMBV,MBSim))
459  std::istringstream iss(strs);
460  char c;
461  iss>>c;
462  if(c=='[') iss.str(strs);
463  else iss.str(std::string("[")+strs+"]");
464 
465  int buf=0;
466  iss >> c;
467  iss >> c;
468  if(c!=']') {
469  iss.putback(c);
470  AT x;
471  do {
472  iss >> x;
473  iss >> c;
474  if(c==';') {
475  if(buf)
476  assert(buf == n);
477 
478  buf=n;
479  n=0;
480  m++;
481  }
482  else if(c==',')
483  n++;
484  c='0';
485  } while(iss);
486 
487  n++; m++;
488  lda=m;
489  tp = false;
490  memory.resize(m*n);
491  ele = (AT*)memory.get();
492  iss.clear();
493  iss.seekg(0);
494  iss >> c;
495  for(int i=0; i<m; i++)
496  for(int j=0; j<n; j++) {
497  iss >> er(i,j);
498  iss >> c;
499  }
500  }
501  }
502 
503  template <class AT>
505 
506  m=A.m;
507  n=A.n;
508  memory = A.memory;
509  ele = A.ele;
510  lda = A.lda;
511  tp = A.tp;
512 
513  return *this;
514  }
515 
516  template <class AT>
518 
519  if(!ele) {
520  m = A.rows();
521  n = A.cols();
522  lda = m;
523  tp = false;
524  memory.resize(m*n);
525  ele = (AT*)memory.get();
526  } else {
527 #ifndef FMATVEC_NO_SIZE_CHECK
528  assert(m == A.rows());
529  assert(n == A.cols());
530 #endif
531  }
532 
533  deepCopy(A);
534 
535  return *this;
536  }
537 
538  template <class AT> template< class Type, class Row, class Col>
540 
541  if(!ele) {
542  m = A.rows();
543  n = A.cols();
544  lda = m;
545  tp = false;
546  memory.resize(m*n);
547  ele = (AT*)memory.get();
548  } else {
549 #ifndef FMATVEC_NO_SIZE_CHECK
550  assert(m == A.rows());
551  assert(n == A.cols());
552 #endif
553  }
554 
555  deepCopy(A);
556 
557  return *this;
558  }
559 
560  template <class AT> template< class Type, class Row, class Col>
562 
563  if(m!=A.rows() || n!=A.cols()) {
564  m = A.rows();
565  n = A.cols();
566  lda = m;
567  tp = false;
568  memory.resize(m*n);
569  ele = (AT*)memory.get();
570  }
571 
572  deepCopy(A);
573 
574  return *this;
575  }
576 
577  template <class AT>
579  for(int i=0; i<m*n; i++) ele[i]=val;
580  return *this;
581  }
582 
583  template <class AT>
585 
586  if(tp) {
587  for(int i=0; i<rows(); i++)
588  for(int j=0; j<cols(); j++)
589  et(i,j) = val;
590  }
591  else {
592  for(int i=0; i<rows(); i++)
593  for(int j=0; j<cols(); j++)
594  er(i,j) = val;
595  }
596 
597  return *this;
598  }
599 
600  template <class AT>
602  if(tp)
603  for(int i=0; i<m; i++)
604  for(int j=0; j<n; j++)
605  et(i,j) = (i==j) ? val : 0;
606  else
607  for(int i=0; i<m; i++)
608  for(int j=0; j<n; j++)
609  er(i,j) = (i==j) ? val : 0;
610  return *this;
611  }
612 
613  template <class AT>
615  return operator()(Index(i1,i2),Index(j1,j2));
616  }
617 
618  template <class AT>
619  inline const Matrix<General,Ref,Ref,AT> Matrix<General,Ref,Ref,AT>::operator()(int i1, int j1, int i2, int j2) const {
620  return operator()(Index(i1,i2),Index(j1,j2));
621  }
622 
623  template <class AT>
625 #ifndef FMATVEC_NO_BOUNDS_CHECK
626  assert(I.end()<m);
627  assert(J.end()<n);
628 #endif
629  return Matrix<General,Ref,Ref,AT>(I.end()-I.start()+1,J.end()-J.start()+1,lda,tp,memory,elePtr(I.start(),J.start()));
630  }
631 
632  template <class AT>
634 #ifndef FMATVEC_NO_BOUNDS_CHECK
635  assert(I.end()<m);
636  assert(J.end()<n);
637 #endif
638  return Matrix<General,Ref,Ref,AT>(I.end()-I.start()+1,J.end()-J.start()+1,lda,tp,memory,elePtr(I.start(),J.start()));
639  }
640 
641  template <class AT>
643 #ifndef FMATVEC_NO_BOUNDS_CHECK
644  assert(I.end()<m);
645 #endif
646  return SquareMatrix<Ref,AT>(I.end()-I.start()+1,lda,tp,memory,elePtr(I.start(),I.start()));
647  }
648 
649  template <class AT>
651 #ifndef FMATVEC_NO_BOUNDS_CHECK
652  assert(I.end()<m);
653 #endif
654  return SquareMatrix<Ref,AT>(I.end()-I.start()+1,lda,tp,memory,elePtr(I.start(),I.start()));
655  }
656 
657  template <class AT>
659 
660 #ifndef FMATVEC_NO_BOUNDS_CHECK
661  assert(i>=0);
662  assert(i<n);
663 #endif
664 
665  return Vector<Ref,AT>(m,lda,tp,memory,elePtr(0,i));
666  }
667 
668  template <class AT>
670 
671 #ifndef FMATVEC_NO_BOUNDS_CHECK
672  assert(i>=0);
673  assert(i<n);
674 #endif
675 
676  return Vector<Ref,AT>(m,lda,tp,memory,elePtr(0,i));
677  }
678 
679  template <class AT>
681 
682 #ifndef FMATVEC_NO_BOUNDS_CHECK
683  assert(i>=0);
684  assert(i<m);
685 #endif
686 
687  return RowVector<Ref,AT>(n,lda,tp,memory,elePtr(i,0));
688  }
689 
690  template <class AT>
692 
693 #ifndef FMATVEC_NO_BOUNDS_CHECK
694  assert(i>=0);
695  assert(i<m);
696 #endif
697 
698  return RowVector<Ref,AT>(n,lda,tp,memory,elePtr(i,0));
699  }
700 
701  template <class AT>
703 
704  Matrix<General,Ref,Ref,AT> A(m,n,NONINIT);
705  A.deepCopy(*this);
706 
707  return A;
708  }
709 
710  template <class AT>
711  inline Matrix<General,Ref,Ref,AT>::operator std::vector<std::vector<AT> >() {
712  std::vector<std::vector<AT> > ret(rows());
713  if(tp) {
714  for(int r=0; r<rows(); r++) {
715  ret[r].resize(cols());
716  for(int c=0; c<cols(); c++)
717  ret[r][c]= et(r,c);
718  }
719  }
720  else {
721  for(int r=0; r<rows(); r++) {
722  ret[r].resize(cols());
723  for(int c=0; c<cols(); c++)
724  ret[r][c]= er(r,c);
725  }
726  }
727  return ret;
728  }
729 
730  template <class AT>
731  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) {
732 #ifndef FMATVEC_NO_INITIALIZATION
733  init(0);
734 #endif
735  for(int r=0; r<rows(); r++) {
736  assert(m[r].size()==cols());
737  for(int c=0; c<cols(); c++)
738  er(r,c)=m[r][c];
739  }
740  }
741 
743 
744  template <class AT> template <class Type, class Row, class Col>
746  if(tp) {
747  for(int i=0; i<m; i++)
748  for(int j=0; j<n; j++)
749  et(i,j) = A.e(i,j);
750  }
751  else {
752  for(int i=0; i<m; i++)
753  for(int j=0; j<n; j++)
754  er(i,j) = A.e(i,j);
755  }
756  }
757 
758  template <class AT>
759  inline void Matrix<General,Ref,Ref,AT>::deepCopy(const Matrix<General,Ref,Ref,AT> &A) {
760  if(A.tp) {
761  if(tp) {
762  for(int i=0; i<m; i++)
763  for(int j=0; j<n; j++)
764  et(i,j) = A.et(i,j);
765  }
766  else {
767  for(int i=0; i<m; i++)
768  for(int j=0; j<n; j++)
769  er(i,j) = A.et(i,j);
770  }
771  } else {
772  if(tp) {
773  for(int i=0; i<m; i++)
774  for(int j=0; j<n; j++)
775  et(i,j) = A.er(i,j);
776  }
777  else {
778  for(int i=0; i<m; i++)
779  for(int j=0; j<n; j++)
780  er(i,j) = A.er(i,j);
781  }
782  }
783  }
784 
785  template<class AT>
786  inline void Matrix<General,Ref,Ref,AT>::deepCopy(const Matrix<Symmetric,Ref,Ref,AT> &A) {
787  for(int i=0; i<A.size(); i++) {
788  er(i,i) = A.ej(i,i);
789  for(int j=i+1; j<A.size(); j++)
790  er(i,j) = et(i,j) = A.ej(i,j);
791  }
792  }
793 
795 
796 }
797 
798 #endif
This is the basic matrix class for arbitrary matrices.
Definition: matrix.h:56
int ldim() const
Leading dimension.
Definition: general_matrix.h:282
AT & operator()(int i, int j)
Standard constructor.
Definition: matrix.h:85
int rows() const
Number of rows.
Definition: general_matrix.h:270
Definition: matrix.h:39
AT & operator()(int i, int j)
Element operator.
Definition: general_matrix.h:203
Matrix(const Matrix< General, Ref, Ref, AT > &A)
Copy Constructor.
Definition: general_matrix.h:104
Definition: matrix.h:218
const CBLAS_TRANSPOSE blasTrans() const
Transposed status.
Definition: general_matrix.h:297
int rows() const
Number of rows.
int end() const
Last element.
Definition: index.h:85
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: matrix.h:212
Definition: matrix.h:215
const CBLAS_ORDER blasOrder() const
Storage convention.
Definition: general_matrix.h:308
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:218
std::istream & operator>>(std::istream &is, Matrix< Type, Row, Col, AT > &A)
Matrix input.
Definition: matrix.h:170
AT * operator()()
Pointer operator.
Definition: general_matrix.h:258
This is an index class for creating submatrices.
Definition: index.h:34
Matrix(int m_, int n_, AT *ele_)
Regular Constructor.
Definition: general_matrix.h:114
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:289
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:38
Matrix()
Standard constructor.
Definition: general_matrix.h:85
int start() const
First element.
Definition: index.h:79
~Matrix()
Destructor.
Definition: general_matrix.h:133
const AT * operator()() const
Pointer operator.
Definition: general_matrix.h:264
int cols() const
Number of columns.
Definition: general_matrix.h:276
Definition: matrix.h:40

Impressum / Disclaimer / Datenschutz Generated by doxygen 1.8.5 Valid HTML