// $Id: tmatrix.cc,v 1.4 2000-04-05 15:44:15 ansari Exp $ // C.Magneville 04/99 #include "machdefs.h" #include #include #include "pexceptions.h" #include "tmatrix.h" //////////////////////////////////////////////////////////////// //**** Createur, Destructeur template TMatrix::TMatrix() // Constructeur par defaut. : TArray() { ck_memo_vt_ = true; } template TMatrix::TMatrix(uint_4 r,uint_4 c, short mm) // Construit une matrice de r lignes et c colonnes. : TArray() { if ( (r == 0) || (c == 0) ) throw ParmError("TMatrix::TMatrix(uint_4 r,uint_4 c) NRows or NCols = 0"); ReSize(r, c, mm); } template TMatrix::TMatrix(const TMatrix& a) // Constructeur par copie (partage si "a" temporaire). : TArray(a) { } template TMatrix::TMatrix(const TMatrix& a, bool share) // Constructeur par copie avec possibilite de forcer le partage ou non. : TArray(a, share) { } template TMatrix::TMatrix(const TArray& a) : TArray(a) { if (a.NbDimensions() > 2) throw SzMismatchError("TMatrix::TMatrix(const TArray& a) a.NbDimensions()>2 "); if (a.NbDimensions() == 1) { size_[1] = 1; step_[1] = size_[0]*step_[0]; ndim_ = 2; } UpdateMemoryMapping(a, SameMemoryMapping); } template TMatrix::TMatrix(const TArray& a, bool share, short mm ) : TArray(a, share) { if (a.NbDimensions() > 2) throw SzMismatchError("TMatrix::TMatrix(const TArray& a, ...) a.NbDimensions()>2"); if (a.NbDimensions() == 1) { size_[1] = 1; step_[1] = size_[0]*step_[0]; ndim_ = 2; } UpdateMemoryMapping(a, mm); } template TMatrix::~TMatrix() // Destructeur { } template TArray& TMatrix::Set(const TArray& a) { if (a.NbDimensions() > 2) throw SzMismatchError("TMatrix::Set(const TArray& a) a.NbDimensions() > 2"); TArray::Set(a); if (a.NbDimensions() == 1) { size_[1] = 1; step_[1] = size_[0]*step_[0]; ndim_ = 2; } UpdateMemoryMapping(a, SameMemoryMapping); return(*this); } template void TMatrix::ReSize(uint_4 r, uint_4 c, short mm) { if(r==0||c==0) throw(SzMismatchError("TMatrix::ReSize r or c==0 ")); uint_4 size[BASEARRAY_MAXNDIMS]; for(int kk=0; kk::ReSize(2, size, 1); UpdateMemoryMapping(mm); } template void TMatrix::Realloc(uint_4 r,uint_4 c, short mm, bool force) { if(r==0||c==0) throw(SzMismatchError("TMatrix::Realloc r or c==0 ")); uint_4 size[BASEARRAY_MAXNDIMS]; for(int kk=0; kk::Realloc(2, size, 1, force); UpdateMemoryMapping(mm); } // $CHECK$ Reza 03/2000 Doit-on declarer cette methode const ? template TMatrix TMatrix::SubMatrix(Range rline, Range rcol) const { short mm = GetMemoryMapping(); Range rx, ry; if (mm == CMemoryMapping) { rx = rcol; ry = rline; } else { ry = rcol; rx = rline; } TMatrix sm(SubArray(rx, ry, Range(0), Range(0), Range(0)),true, mm); sm.UpdateMemoryMapping(mm); sm.SetTemp(true); return(sm); } //////////////////////////////////////////////////////////////// // Transposition template TMatrix& TMatrix::Transpose() { short vt = (marowi_ == veceli_) ? ColumnVector : RowVector; uint_4 rci = macoli_; macoli_ = marowi_; marowi_ = rci; veceli_ = (vt == ColumnVector ) ? marowi_ : macoli_; return(*this); } template TMatrix TMatrix::Transpose(short mm) { if (mm == SameMemoryMapping) mm = GetMemoryMapping(); TMatrix tm(NCols(), NRows(), mm); for(uint_4 i=0; i TMatrix TMatrix::Rearrange(short mm) { if ( mm == SameMemoryMapping) mm = GetMemoryMapping(); else if ( (mm != CMemoryMapping) && (mm != FortranMemoryMapping) ) mm = GetDefaultMemoryMapping(); if (mm == GetMemoryMapping()) return (TMatrix(*this, true)); TMatrix tm(NRows(), NCols(), mm); for(uint_4 i=0; i TMatrix& TMatrix::SetIdentity(IdentityMatrix imx) { if (ndim_ == 0) { uint_4 sz = imx.Size(); if (sz < 1) sz = 1; ReSize(sz, sz); } T diag = (T)imx.Diag(); if (NRows() != NCols()) throw SzMismatchError("TMatrix::operator= (IdentityMatrix) NRows() != NCols()") ; for(uint_4 i=0; i string TMatrix::InfoString() const { string rs = "TMatrix<"; rs += typeid(T).name(); char buff[64]; sprintf(buff, ">(NRows=%ld, NCols=%ld)", (long)NRows(), (long)NCols()); rs += buff; return(rs); } template void TMatrix::Print(ostream& os, int_4 maxprt, bool si) const { if (maxprt < 0) maxprt = max_nprt_; int_4 npr = 0; Show(os, si); uint_4 kc,kr; for(kr=0; kr 1) && (size_[macoli_] > 10) ) cout << "----- Ligne Line= " << kr << endl; for(kc=0; kc 0) os << ", "; os << (*this)(kr, kc); npr++; if (npr >= maxprt) { if (npr < totsize_) os << "\n .... " << endl; return; } } os << endl; } os << endl; } //////////////////////////////////////////////////////////////// //**** Multiplication matricielle ***** //////////////////////////////////////////////////////////////// template TMatrix TMatrix::Multiply(const TMatrix& b, short mm) const { if (NCols() != b.NRows()) throw(SzMismatchError("TMatrix::Multiply(b) NCols() != b.NRows() ") ); if (mm == SameMemoryMapping) mm = GetMemoryMapping(); TMatrix rm(NRows(), b.NCols(), mm); const T * pea; const T * peb; T sum; uint_4 r,c,k; uint_4 stepa = Step(ColsKA()); uint_4 stepb = b.Step(RowsKA()); // Calcul de C=rm = A*B (A=*this) for(r=0; r #pragma define_template TMatrix #pragma define_template TMatrix #pragma define_template TMatrix #pragma define_template TMatrix #pragma define_template TMatrix< complex > #pragma define_template TMatrix< complex > #endif #if defined(ANSI_TEMPLATES) || defined(GNU_TEMPLATES) template class TMatrix; template class TMatrix; template class TMatrix; template class TMatrix; template class TMatrix; template class TMatrix< complex >; template class TMatrix< complex >; #endif