Changeset 286 in Sophya


Ignore:
Timestamp:
Apr 30, 1999, 1:02:52 PM (26 years ago)
Author:
ansari
Message:

tmatrix suite et fin cmv 30/4/99

Location:
trunk/SophyaLib/NTools
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • trunk/SophyaLib/NTools/tmatrix.cc

    r279 r286  
    1 // $Id: tmatrix.cc,v 1.1 1999-04-29 08:57:15 ansari Exp $
     1// $Id: tmatrix.cc,v 1.2 1999-04-30 11:02:52 ansari Exp $
    22//                         C.Magneville          04/99
    33#include "machdefs.h"
    … …  
    1818TMatrix<T>::TMatrix()
    1919// Constructeur par defaut.
    20 : mNr(0), mNc(0), mNDBlock(NULL)
     20: mNr(0), mNc(0), mNDBlock()
    2121{
    2222}
    … …  
    2525TMatrix<T>::TMatrix(uint_4 r,uint_4 c)
    2626// Construit une matrice de r lignes et c colonnes.
    27 : mNr(r), mNc(c), mNDBlock(new NDataBlock<T>(r*c))
     27: mNr(r), mNc(c), mNDBlock(r*c)
    2828{
    2929}
    … …  
    3333// Construit une matrice de r lignes et c colonnes. On fournit
    3434// le tableau des valeurs et eventuellement un Bridge.
    35 : mNr(r), mNc(c), mNDBlock(new NDataBlock<T>(r*c,values,br))
     35: mNr(r), mNc(c), mNDBlock(r*c,values,br)
    3636{
    3737}
    … …  
    3939template <class T>
    4040TMatrix<T>::TMatrix(const TMatrix<T>& a)
    41 // Constructeur par copie.
    42 : mNr(a.mNr), mNc(a.mNc)
    43 , mNDBlock(new NDataBlock<T>(*(a.DataBlock())))
     41// Constructeur par copie (partage des donnees).
     42: mNr(a.mNr), mNc(a.mNc), mNDBlock(a.mNDBlock)
    4443{
    4544}
    … …  
    4847TMatrix<T>::TMatrix(const TMatrix<T>& a,bool share)
    4948// Constructeur par copie.
    50 : mNr(a.mNr), mNc(a.mNc)
    51 , mNDBlock(new NDataBlock<T>(*(a.DataBlock()),share))
     49: mNr(a.mNr), mNc(a.mNc), mNDBlock(a.mNDBlock,share)
    5250{
    5351}
    … …  
    5755// Destructeur
    5856{
    59 Delete();
    60 }
    61 
    62 template <class T>
    63 void TMatrix<T>::Delete(void)
    64 // Pour desallouer
    65 {
    66 mNr = mNc = 0;
    67 if(mNDBlock) delete mNDBlock; mNDBlock = NULL;
    6857}
    6958
    … …  
    7160template <class T>
    7261void TMatrix<T>::Clone(const TMatrix<T>& a)
    73 // Clone (copie de donnee) a partir de "a"
    74 {
    75 if(a.mNr==0 || a.mNc==0 || a.DataBlock()==NULL) {Delete(); return;}
    76 if(mNDBlock) Delete();
    77 else NDataBlock<T>* mNDBlock = new NDataBlock<T>(1); // cas "a" vide
    78 mNr = a.mNr; mNc = a.mNc;
    79 mNDBlock->Clone(*(a.DataBlock()));
     62// Clone (copie de donnees) a partir de "a"
     63// Partage des donnees si "a" est temporaire.
     64{
     65mNDBlock.Clone(a.mNDBlock); mNr = a.mNr; mNc = a.mNc;
    8066}
    8167
    … …  
    8470// Reset de la matrice a "v"
    8571{
    86 if(mNDBlock==NULL)
    87   throw(NullPtrError("TMatrix::Reset mNDBlock==NULL\n"));
    88 if(mNr==0 || mNc==0)
    89   throw(SzMismatchError("TMatrix::Reset mNr==0 || mNc==0\n"));
    90 mNDBlock->Reset(v);
     72mNDBlock.Reset(v);
    9173}
    9274
    … …  
    9678{
    9779if(r==0 || c==0) throw(SzMismatchError("TMatrix::ReSize r ou c==0\n"));
    98 if(mNDBlock==NULL) NDataBlock<T>* mNDBlock = new NDataBlock<T>(1);
    9980if(!force_alloc && mNr==r && mNc==c)  return;
    100 mNr = r; mNc = c;
    101 mNDBlock->ReSize(r*c);
     81mNr = r; mNc = c; mNDBlock.ReSize(r*c,force_alloc);
    10282}
    10383
    … …  
    10989// Operateur d'affectation depuis scalaire : identite * scalaire.
    11090{
    111 if(mNr!=mNc || mNr==0 || mNDBlock==NULL)
    112   throw(SzMismatchError("TMatrix::operator= mNc!=mNr ou ==0 ou mNDBlock==NULL\n"));
    113 for(int r=0;r<mNr;r++) for(int c=0;c<mNc;c++) (*this)(r,c) = (r==c)? x: 0;
     91if(mNr!=mNc || mNr==0)
     92  throw(SzMismatchError("TMatrix::operator= mNc!=mNr ou ==0\n"));
     93for(uint_4 r=0;r<mNr;r++) for(uint_4 c=0;c<mNc;c++)
     94                               (*this)(r,c) = (r==c)? x: 0;
    11495return *this;
    11596}
    … …  
    11798template <class T>
    11899TMatrix<T>& TMatrix<T>::operator = (const TMatrix<T>& a)
    119 // Operateur d'affectation A=B (sans reaffectation)
     100// Operateur d'affectation A=a
     101// (Re)allocation de la place memoire
     102// sauf si "a" est cree temporairement (partage)
    120103{
    121104if(this == &a) return *this;
    122 if(mNr!=a.mNc || mNc!=a.mNc || mNr==0 || mNc==0)
    123   throw(SzMismatchError("TMatrix::operator= mNc!=a.mNc ou mNr!=a.mNr ou ==0\n"));
    124 if(mNDBlock==NULL)
    125   throw(SzMismatchError("TMatrix::operator= mNDBlock==NULL\n"));
    126 *(mNDBlock) = *(a.DataBlock());
     105Clone(a);
    127106return *this;
    128107}
    … …  
    138117  os<<"TMatrix::Print("<<mNr<<","<<mNc<<")"<<endl;
    139118if(lp>0)
    140   {os<<"       this="<<this<<" mNDBlock="<<mNDBlock<<endl;
    141    if(mNDBlock) mNDBlock->Print(0,0);}
     119  {os<<"       this="<<this<<endl; mNDBlock.Print(0,0);}
    142120if(mNr==0 || mNc==0) return;
    143121if(i0>=mNr || j0>=mNc || ni==0 || nj==0) return;
    144122uint_4 i1 = i0+ni; if(i1>mNr) i1 = mNr;
    145123uint_4 j1 = j0+nj; if(j1>mNc) j1 = mNc;
    146 for(int i=i0;i<i1;i++) {
    147   for(int j=j0;j<j1;j++) cout<<" "<<(*this)(i,j)<<endl;
     124for(uint_4 i=i0;i<i1;i++) {
     125  for(uint_4 j=j0;j<j1;j++) cout<<" "<<(*this)(i,j);
     126  cout<<endl;
    148127}
    149128}
    … …  
    155134TMatrix<T>& TMatrix<T>::operator += (T b)
    156135{
    157 if(mNr==0 || mNc==0 || !mNDBlock)
    158   throw(SzMismatchError("TMatrix::operator+=v mNr==0 || mNc==0 || !mNDBlock\n"));
    159 *mNDBlock += b;
     136mNDBlock += b;
    160137return *this;
    161138}
    … …  
    164141TMatrix<T>& TMatrix<T>::operator -= (T b)
    165142{
    166 if(mNr==0 || mNc==0 || !mNDBlock)
    167   throw(SzMismatchError("TMatrix::operator-=v mNr==0 || mNc==0 || !mNDBlock\n"));
    168 *mNDBlock -= b;
     143mNDBlock -= b;
    169144return *this;
    170145}
    … …  
    173148TMatrix<T>& TMatrix<T>::operator *= (T b)
    174149{
    175 if(mNr==0 || mNc==0 || !mNDBlock)
    176   throw(SzMismatchError("TMatrix::operator*=v mNr==0 || mNc==0 || !mNDBlock\n"));
    177 *mNDBlock *= b;
     150mNDBlock *= b;
    178151return *this;
    179152}
    … …  
    182155TMatrix<T>& TMatrix<T>::operator /= (T b)
    183156{
    184 if(b==(T) 0) throw(ParmError("TMatrix::operator/=v divide by zero\n"));
    185 if(mNr==0 || mNc==0 || !mNDBlock)
    186   throw(SzMismatchError("TMatrix::operator/=v mNr==0 || mNc==0 || !mNDBlock\n"));
    187 *mNDBlock /= b;
     157mNDBlock /= b;
    188158return *this;
    189159}
    … …  
    195165TMatrix<T>& TMatrix<T>::operator += (const TMatrix<T>& a)
    196166{
    197 if(mNr==0 || mNc==0 || mNr!=a.mNr || mNc!=a.mNc || !mNDBlock || !a.mNDBlock)
    198   throw(SzMismatchError("TMatrix::operator+=A size mismatch/null"));
    199 *(mNDBlock) += *(a.mNDBlock);
     167if(mNr==0 || mNc==0 || mNr!=a.mNr || mNc!=a.mNc)
     168  throw(SzMismatchError("TMatrix::operator+=A size mismatch"));
     169mNDBlock += a.mNDBlock;
    200170return *this;
    201171}
    … …  
    204174TMatrix<T>& TMatrix<T>::operator -= (const TMatrix<T>& a)
    205175{
    206 if(mNr==0 || mNc==0 || mNr!=a.mNr || mNc!=a.mNc || !mNDBlock || !a.mNDBlock)
    207   throw(SzMismatchError("TMatrix::operator-=A size mismatch/null"));
    208 *(mNDBlock) -= *(a.mNDBlock);
     176if(mNr==0 || mNc==0 || mNr!=a.mNr || mNc!=a.mNc)
     177  throw(SzMismatchError("TMatrix::operator-=A size mismatch"));
     178mNDBlock -= a.mNDBlock;
    209179return *this;
    210180}
    … …  
    215185{
    216186uint_4 ndata = mNr*mNc;
    217 if(ndata==0 || mNc!=a.mNr || a.mNr!=a.mNc || !mNDBlock || !a.mNDBlock)
    218   throw(SzMismatchError("TMatrix::operator*=A size mismatch/null"));
     187if(ndata==0 || mNc!=a.mNr || a.mNr!=a.mNc)
     188  throw(SzMismatchError("TMatrix::operator*=A size mismatch"));
    219189// A(i,j) = Sum(k) A(i,k)*B(k,j) ... il faut sauver la ligne "i" de A
    220190// Vecteur oi   : vecteur ou est sauve la ligne i de la matrice *this
    … …  
    238208}
    239209
     210////////////////////////////////////////////////////////////////
     211//**** Surcharge de +,-,*: TMatrix = TMatrix + TMatrix;
     212
     213template <class T>
     214TMatrix<T> TMatrix<T>::Add(const TMatrix<T>& b) const
     215// C = A(this)+B
     216{
     217if(mNr!=b.mNr || mNc!=b.mNc)
     218  throw(SzMismatchError("NDataBlock operator C=A+B size mismatch\n"));
     219TMatrix<T> result; result.SetTemp(true);
     220if(b.IsTemp()) {
     221  result.Clone(b);
     222  result += *this;
     223} else {
     224  result.Clone(*this);
     225  result += b;
     226}
     227return result;
     228}
     229
     230template <class T>
     231TMatrix<T> TMatrix<T>::Sub(const TMatrix<T>& b) const
     232// C = A(this)-B
     233{
     234if(mNr!=b.mNr || mNc!=b.mNc)
     235  throw(SzMismatchError("NDataBlock operator C=A-B size mismatch\n"));
     236TMatrix<T> result; result.SetTemp(true);
     237if(b.IsTemp()) {
     238  result.Clone(b);
     239  result.mNDBlock = (*this).mNDBlock - result.mNDBlock;
     240} else {
     241  result.Clone(*this);
     242  result -= b;
     243}
     244return result;
     245}
     246
     247template <class T>
     248TMatrix<T> TMatrix<T>::Mul(const TMatrix<T>& b) const
     249// C = A(this)*B : Cij = Aik Bkj (allocation forcee dans tous les cas)
     250{
     251if(mNr==0 || mNc==0 || b.mNc==0 || mNc!=b.mNr)
     252  throw(SzMismatchError("NDataBlock operator C=A*B size mismatch\n"));
     253TMatrix<T> r; r.SetTemp(true); r.Realloc(mNr,b.mNc,true);
     254T *ai,*aik,*bj,*bkj,*ri,*rij;
     255for(ri=const_cast<T *>(r.Data()),ai=const_cast<T *>(Data());
     256    ri<r.Data()+r.mNr*r.mNc;ri+=r.mNc,ai+=mNc) {
     257  for(rij=ri,bj=const_cast<T *>(b.Data());rij<ri+r.mNc;rij++,bj++) {
     258    *rij = 0;
     259    for(aik=ai,bkj=bj;aik<ai+mNc;aik++,bkj+=b.mNc) *rij += *aik * *bkj;
     260  }
     261}
     262return r;
     263}
     264
     265////////////////////////////////////////////////////////////////
     266// -------------------------------------------------------------------------
     267//   Les objets delegues pour la gestion de persistance
     268// -------------------------------------------------------------------------
     269
     270template <class T>
     271FIO_TMatrix<T>::FIO_TMatrix()
     272{
     273dobj=new TMatrix<T>;
     274ownobj=true;
     275}
     276
     277template <class T>
     278FIO_TMatrix<T>::FIO_TMatrix(string const & filename)
     279{
     280dobj=new TMatrix<T>;
     281ownobj=true;
     282Read(filename);
     283}
     284
     285template <class T>
     286FIO_TMatrix<T>::FIO_TMatrix(const TMatrix<T> & obj)
     287{
     288dobj = new TMatrix<T>(obj);
     289ownobj=true;
     290}
     291
     292template <class T>
     293FIO_TMatrix<T>::FIO_TMatrix(TMatrix<T> * obj)
     294{
     295dobj = obj;
     296ownobj=false;
     297}
     298
     299template <class T>
     300FIO_TMatrix<T>::~FIO_TMatrix()
     301{
     302if (ownobj && dobj) delete dobj;
     303}
     304
     305template <class T>
     306AnyDataObj* FIO_TMatrix<T>::DataObj()
     307{
     308return(dobj);
     309}
     310
     311
     312template <class T>
     313void FIO_TMatrix<T>::ReadSelf(PInPersist& is)
     314{
     315// On lit les 3 premiers uint_8
     316uint_4 itab[2];
     317is.Get(itab,2);
     318if (dobj == NULL) dobj = new TMatrix<T>(itab[0],itab[1]);
     319else dobj->Realloc(itab[0],itab[1],false);
     320// On lit le NDataBlock
     321//cmv   encore des problemes avec les consteries
     322//FIO_NDataBlock<T> fio_nd(&dobj->DataBlock());
     323//fio_nd.ReadSelf(is);
     324}
     325
     326template <class T>
     327void FIO_TMatrix<T>::WriteSelf(POutPersist& os) const
     328{
     329if (dobj == NULL)   return;
     330//  On ecrit 2 uint_4 .... 0 : NRows,  1 : NCol
     331uint_4 itab[2];
     332itab[0] = dobj->NRows();
     333itab[1] = dobj->NCol();
     334os.Put(itab,2);
     335// On ecrit le NDataBlock
     336//cmv   encore des problemes avec les consteries
     337//FIO_NDataBlock<T> fio_nd(&dobj->DataBlock());
     338//fio_nd.WriteSelf(os);
     339}
     340
    240341///////////////////////////////////////////////////////////////
    241342#ifdef __CXX_PRAGMA_TEMPLATES__
    242 /*
    243343#pragma define_template TMatrix<uint_1>
    244344#pragma define_template TMatrix<uint_2>
    … …  
    249349#pragma define_template TMatrix<uint_8>
    250350#pragma define_template TMatrix<r_4>
    251 */
    252351#pragma define_template TMatrix<r_8>
    253 /*
    254 #pragma define_template TMatrix< complex<float> >
    255 #pragma define_template TMatrix< complex<double> >
    256 */
     352#pragma define_template TMatrix< complex<float> >
     353#pragma define_template TMatrix< complex<double> >
    257354// Instances des delegues FileIO (PPersist)
     355#pragma define_template FIO_TMatrix<uint_1>
     356#pragma define_template FIO_TMatrix<uint_2>
     357#pragma define_template FIO_TMatrix<int_2>
     358#pragma define_template FIO_TMatrix<int_4>
     359#pragma define_template FIO_TMatrix<int_8>
     360#pragma define_template FIO_TMatrix<uint_4>
     361#pragma define_template FIO_TMatrix<uint_8>
     362#pragma define_template FIO_TMatrix<r_8>
     363#pragma define_template FIO_TMatrix<r_4>
     364#pragma define_template FIO_TMatrix< complex<float> >
     365#pragma define_template FIO_TMatrix< complex<double> >
    258366#endif
    259 
    260367
    261368#if defined(ANSI_TEMPLATES) || defined(GNU_TEMPLATES)
    … …  
    272379template class TMatrix< complex<double> >;
    273380// Instances des delegues FileIO (PPersist)
     381template class FIO_TMatrix<uint_1>;
     382template class FIO_TMatrix<uint_2>;
     383template class FIO_TMatrix<int_2>;
     384template class FIO_TMatrix<int_4>;
     385template class FIO_TMatrix<int_8>;
     386template class FIO_TMatrix<uint_4>;
     387template class FIO_TMatrix<uint_8>;
     388template class FIO_TMatrix<r_8>;
     389template class FIO_TMatrix<r_4>;
     390template class FIO_TMatrix< complex<float> >;
     391template class FIO_TMatrix< complex<double> >;
    274392#endif

  • trunk/SophyaLib/NTools/tmatrix.h

    r279 r286  
    2424  virtual ~TMatrix();
    2525
     26  // Temporaire?
     27  inline bool IsTemp(void) const {return mNDBlock.IsTemp();}
     28  inline void SetTemp(bool temp=false) const
     29                     {mNDBlock.SetTemp(temp);}
     30
    2631  // Gestion taille/Remplissage
    2732  void Clone(const TMatrix<T>& a);
    … …  
    3338  inline int NCol()  const {return mNc;}
    3439  T const& operator()(uint_4 r,uint_4 c) const
    35                      {return *(mNDBlock->Begin()+r*mNc+c);}
     40                     {return *(mNDBlock.Begin()+r*mNc+c);}
    3641  T&       operator()(uint_4 r,uint_4 c)
    37                      {return *(mNDBlock->Begin()+r*mNc+c);}
    38   inline       T* Data()       {return mNDBlock->Begin();}
    39   inline const T* Data() const {return mNDBlock->Begin();}
    40   inline       NDataBlock<T>* DataBlock()       {return mNDBlock;}
    41   inline const NDataBlock<T>* DataBlock() const {return mNDBlock;}
     42                     {return *(mNDBlock.Begin()+r*mNc+c);}
     43  inline       T* Data()       {return mNDBlock.Begin();}
     44  inline const T* Data() const {return mNDBlock.Begin();}
     45  inline       NDataBlock<T>& DataBlock()       {return mNDBlock;}
     46  inline const NDataBlock<T>& DataBlock() const {return mNDBlock;}
    4247
    4348  // Operateur d'affectation
    … …  
    5055              {Print(cout,lp,i0,ni,j0,nj);}
    5156
    52   // Surcharge d'operateurs
     57  // Surcharge d'operateurs INPLACE
    5358  TMatrix<T>& operator += (T b);
    5459  TMatrix<T>& operator -= (T b);
    … …  
    6065  TMatrix<T>& operator *= (const TMatrix<T>& a);
    6166
     67  // Pour surcharge d'operateurs
     68  TMatrix<T> Add(const TMatrix<T>& b) const;
     69  TMatrix<T> Sub(const TMatrix<T>& b) const;
     70  TMatrix<T> Mul(const TMatrix<T>& b) const;
     71
    6272protected:
    63   void Delete(void);
    64 
    6573  uint_4 mNr,mNc;
    66   NDataBlock<T>* mNDBlock;
     74  NDataBlock<T> mNDBlock;
    6775};
    6876
     77/////////////////////////////////////////////////////////////////////////
    6978template<class T>
    7079inline ostream& operator << (ostream& os, const TMatrix<T>& a)
    7180                      {a.Print(os); return(os);}
    7281
     82/////////////////////////////////////////////////////////////////////////
     83template<class T>
     84inline TMatrix<T> operator + (const TMatrix<T>& a, double b) {
     85  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
     86  result += b;
     87  return result;
     88}
     89
     90template<class T>
     91inline TMatrix<T> operator + (double b,const TMatrix<T>& a) {
     92  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
     93  result += b;
     94  return result;
     95}
     96
     97template<class T>
     98inline TMatrix<T> operator - (const TMatrix<T>& a, double b) {
     99  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
     100  result -= b;
     101  return result;
     102}
     103
     104template<class T>
     105inline TMatrix<T> operator - (double b,const TMatrix<T>& a) {
     106  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
     107  result.DataBlock() = b - result.DataBlock();
     108  return result;
     109}
     110
     111template<class T>
     112inline TMatrix<T> operator * (const TMatrix<T>& a, double b) {
     113  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
     114  result *= b;
     115  return result;
     116}
     117
     118template<class T>
     119inline TMatrix<T> operator * (double b,const TMatrix<T>& a) {
     120  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
     121  result *= b;
     122  return result;
     123}
     124
     125template<class T>
     126inline TMatrix<T> operator / (const TMatrix<T>& a, double b) {
     127  TMatrix<T> result(a,false); result.DataBlock().SetTemp(true);
     128  result /= b;
     129  return result;
     130}
     131
     132/////////////////////////////////////////////////////////////////////////
     133template<class T>
     134inline TMatrix<T> operator + (const TMatrix<T>& a,const TMatrix<T>& b)
     135                             {return a.Add(b);}
     136
     137template<class T>
     138inline TMatrix<T> operator - (const TMatrix<T>& a,const TMatrix<T>& b)
     139                             {return a.Sub(b);}
     140
     141template<class T>
     142inline TMatrix<T> operator * (const TMatrix<T>& a,const TMatrix<T>& b)
     143                             {return a.Mul(b);}
     144
     145/////////////////////////////////////////////////////////////////////////
     146// Classe pour la gestion de persistance
     147template <class T>
     148class FIO_TMatrix : public  PPersist  {
     149
     150public:
     151            FIO_TMatrix();
     152            FIO_TMatrix(string const & filename);
     153            FIO_TMatrix(const TMatrix<T> & obj);
     154            FIO_TMatrix(TMatrix<T> * obj);
     155  virtual   ~FIO_TMatrix();
     156
     157  virtual   AnyDataObj* DataObj();
     158  inline operator TMatrix<T>() { return(*dobj); }
     159
     160protected :
     161  virtual void       ReadSelf(PInPersist&);           
     162  virtual void       WriteSelf(POutPersist&) const; 
     163  TMatrix<T> * dobj;
     164  bool ownobj;
     165};
     166
    73167} // Fin du namespace
    74168
Note: See TracChangeset for help on using the changeset viewer.