Back to home page

EIC code displayed by LXR
 
 

    

Warning, file /include/root/TMatrixTSparse.h was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).

0001 // @(#)root/matrix:$Id$
0002 // Authors: Fons Rademakers, Eddy Offermann   Feb 2004
0003 
0004 /*************************************************************************
0005  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers.               *
0006  * All rights reserved.                                                  *
0007  *                                                                       *
0008  * For the licensing terms see $ROOTSYS/LICENSE.                         *
0009  * For the list of contributors see $ROOTSYS/README/CREDITS.             *
0010  *************************************************************************/
0011 
0012 #ifndef ROOT_TMatrixTSparse
0013 #define ROOT_TMatrixTSparse
0014 
0015 #include "TMatrixTBase.h"
0016 #include "TMatrixTUtils.h"
0017 
0018 #include <cstring>
0019 
0020 #ifdef CBLAS
0021 #include <vecLib/vBLAS.h>
0022 //#include <cblas.h>
0023 #endif
0024 
0025 //////////////////////////////////////////////////////////////////////////
0026 //                                                                      //
0027 // TMatrixTSparse                                                       //
0028 //                                                                      //
0029 // Template class of a general sparse matrix in the Harwell-Boeing      //
0030 // format                                                               //
0031 //                                                                      //
0032 //////////////////////////////////////////////////////////////////////////
0033 
0034 template<class Element> class TMatrixT;
0035 
0036 template<class Element> class TMatrixTSparse : public TMatrixTBase<Element> {
0037 
0038 protected:
0039 
0040    Int_t   *fRowIndex;  //[fNrowIndex] row index
0041    Int_t   *fColIndex;  //[fNelems]    column index
0042    Element *fElements;  //[fNelems]
0043 
0044    void Allocate(Int_t nrows,Int_t ncols,Int_t row_lwb = 0,Int_t col_lwb = 0,
0045                  Int_t init = 0,Int_t nr_nonzeros = 0);
0046 
0047   // Elementary constructors
0048    void AMultB (const TMatrixTSparse<Element> &a,const TMatrixTSparse<Element> &b,Int_t constr=0) {
0049                 const TMatrixTSparse<Element> bt(TMatrixTSparse::kTransposed,b); AMultBt(a,bt,constr); }
0050    void AMultB (const TMatrixTSparse<Element> &a,const TMatrixT<Element>       &b,Int_t constr=0) {
0051                 const TMatrixTSparse<Element> bsp = b;
0052                 const TMatrixTSparse<Element> bt(TMatrixTSparse::kTransposed,bsp); AMultBt(a,bt,constr); }
0053    void AMultB (const TMatrixT<Element>       &a,const TMatrixTSparse<Element> &b,Int_t constr=0) {
0054                 const TMatrixTSparse<Element> bt(TMatrixTSparse::kTransposed,b); AMultBt(a,bt,constr); }
0055 
0056    void AMultBt(const TMatrixTSparse<Element> &a,const TMatrixTSparse<Element> &b,Int_t constr=0);
0057    void AMultBt(const TMatrixTSparse<Element> &a,const TMatrixT<Element>       &b,Int_t constr=0);
0058    void AMultBt(const TMatrixT<Element>       &a,const TMatrixTSparse<Element> &b,Int_t constr=0);
0059 
0060    void APlusB (const TMatrixTSparse<Element> &a,const TMatrixTSparse<Element> &b,Int_t constr=0);
0061    void APlusB (const TMatrixTSparse<Element> &a,const TMatrixT<Element>       &b,Int_t constr=0);
0062    void APlusB (const TMatrixT<Element>       &a,const TMatrixTSparse<Element> &b,Int_t constr=0) { APlusB(b,a,constr); }
0063 
0064    void AMinusB(const TMatrixTSparse<Element> &a,const TMatrixTSparse<Element> &b,Int_t constr=0);
0065    void AMinusB(const TMatrixTSparse<Element> &a,const TMatrixT<Element>       &b,Int_t constr=0);
0066    void AMinusB(const TMatrixT<Element>       &a,const TMatrixTSparse<Element> &b,Int_t constr=0);
0067 
0068 public:
0069 
0070    enum EMatrixCreatorsOp1 { kZero,kUnit,kTransposed,kAtA };
0071    enum EMatrixCreatorsOp2 { kMult,kMultTranspose,kPlus,kMinus };
0072 
0073    TMatrixTSparse() { fElements = nullptr; fRowIndex = nullptr; fColIndex = nullptr; }
0074    TMatrixTSparse(Int_t nrows,Int_t ncols);
0075    TMatrixTSparse(Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb);
0076    TMatrixTSparse(Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb,Int_t nr_nonzeros,
0077                   Int_t *row, Int_t *col,Element *data);
0078    TMatrixTSparse(const TMatrixTSparse<Element> &another);
0079    TMatrixTSparse(const TMatrixT<Element>       &another);
0080 
0081    TMatrixTSparse(EMatrixCreatorsOp1 op,const TMatrixTSparse<Element> &prototype);
0082    TMatrixTSparse(const TMatrixTSparse<Element> &a,EMatrixCreatorsOp2 op,const TMatrixTSparse<Element> &b);
0083    TMatrixTSparse(const TMatrixTSparse<Element> &a,EMatrixCreatorsOp2 op,const TMatrixT      <Element> &b);
0084    TMatrixTSparse(const TMatrixT      <Element> &a,EMatrixCreatorsOp2 op,const TMatrixTSparse<Element> &b);
0085 
0086    ~TMatrixTSparse() override { TMatrixTSparse::Clear(); }
0087 
0088    const Element *GetMatrixArray  () const override;
0089          Element *GetMatrixArray  () override;
0090    const Int_t    *GetRowIndexArray() const override;
0091          Int_t    *GetRowIndexArray() override;
0092    const Int_t    *GetColIndexArray() const override;
0093          Int_t    *GetColIndexArray() override;
0094 
0095    TMatrixTBase<Element>   &SetRowIndexArray(Int_t *data) override { memmove(fRowIndex,data,(this->fNrows+1)*sizeof(Int_t)); return *this; }
0096    TMatrixTBase<Element>   &SetColIndexArray(Int_t *data) override { memmove(fColIndex,data,this->fNelems*sizeof(Int_t)); return *this; }
0097 
0098            TMatrixTSparse<Element> &SetSparseIndex  (Int_t nelem_new);
0099            TMatrixTSparse<Element> &SetSparseIndex  (const TMatrixTBase<Element> &another);
0100            TMatrixTSparse<Element> &SetSparseIndexAB(const TMatrixTSparse<Element> &a,const TMatrixTSparse<Element> &b);
0101            TMatrixTSparse<Element> &SetSparseIndexAB(const TMatrixT      <Element> &a,const TMatrixTSparse<Element> &b);
0102            TMatrixTSparse<Element> &SetSparseIndexAB(const TMatrixTSparse<Element> &a,const TMatrixT      <Element> &b)
0103                                               { return SetSparseIndexAB(b,a); }
0104 
0105    void                     GetMatrix2Array (Element *data,Option_t * /*option*/ ="") const override;
0106    TMatrixTBase<Element>   &SetMatrixArray  (const Element *data,Option_t * /*option*/="") override
0107                                                     { memcpy(fElements,data,this->fNelems*sizeof(Element)); return *this; }
0108    virtual TMatrixTBase<Element>   &SetMatrixArray  (Int_t nr_nonzeros,Int_t *irow,Int_t *icol,Element *data);
0109    virtual TMatrixTBase<Element> &
0110    SetMatrixArray(Int_t nr_nonzeros, Int_t nrows, Int_t ncols, Int_t *irow, Int_t *icol, Element *data);
0111    TMatrixTBase<Element>   &InsertRow       (Int_t row,Int_t col,const Element *v,Int_t n=-1) override;
0112    void                     ExtractRow      (Int_t row,Int_t col,      Element *v,Int_t n=-1) const override;
0113 
0114    TMatrixTBase<Element>   &ResizeTo(Int_t nrows,Int_t ncols,Int_t nr_nonzeros=-1) override;
0115    TMatrixTBase<Element>   &ResizeTo(Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb,Int_t nr_nonzeros=-1) override;
0116    inline  TMatrixTBase<Element>   &ResizeTo(const TMatrixTSparse<Element> &m) {return ResizeTo(m.GetRowLwb(),m.GetRowUpb(),m.GetColLwb(),
0117                                                                                                 m.GetColUpb(),m.GetNoElements()); }
0118 
0119    void Clear(Option_t * /*option*/ ="") override { if (this->fIsOwner) {
0120                                                       if (fElements) { delete [] fElements; fElements = nullptr; }
0121                                                       if (fRowIndex) { delete [] fRowIndex; fRowIndex = nullptr; }
0122                                                       if (fColIndex) { delete [] fColIndex; fColIndex = nullptr; }
0123                                                    }
0124                                                    this->fNelems    = 0;
0125                                                    this->fNrowIndex = 0;
0126                                                  }
0127 
0128            TMatrixTSparse<Element> &Use   (Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb,Int_t nr_nonzeros,
0129                                            Int_t *pRowIndex,Int_t *pColIndex,Element *pData);
0130    const   TMatrixTSparse<Element> &Use   (Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb,Int_t nr_nonzeros,
0131                                            const Int_t *pRowIndex,const Int_t *pColIndex,const Element *pData) const
0132                                             { return (const TMatrixTSparse<Element>&)
0133                                                      ((const_cast<TMatrixTSparse<Element> *>(this))->Use(row_lwb,row_upb,col_lwb,col_upb,nr_nonzeros,
0134                                                                                              const_cast<Int_t *>(pRowIndex),
0135                                                                                              const_cast<Int_t *>(pColIndex),
0136                                                                                              const_cast<Element *>(pData))); }
0137            TMatrixTSparse<Element> &Use   (Int_t nrows,Int_t ncols,Int_t nr_nonzeros,
0138                                            Int_t *pRowIndex,Int_t *pColIndex,Element *pData);
0139    const   TMatrixTSparse<Element> &Use   (Int_t nrows,Int_t ncols,Int_t nr_nonzeros,
0140                                            const Int_t *pRowIndex,const Int_t *pColIndex,const Element *pData) const;
0141            TMatrixTSparse<Element> &Use   (TMatrixTSparse<Element> &a);
0142    const   TMatrixTSparse<Element> &Use   (const TMatrixTSparse<Element> &a) const;
0143 
0144    TMatrixTBase<Element>   &GetSub(Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb,
0145                                             TMatrixTBase<Element> &target,Option_t *option="S") const override;
0146            TMatrixTSparse<Element>  GetSub(Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb,Option_t *option="S") const;
0147    TMatrixTBase<Element>   &SetSub(Int_t row_lwb,Int_t col_lwb,const TMatrixTBase<Element> &source) override;
0148 
0149    Bool_t IsSymmetric() const override { return (*this == TMatrixTSparse<Element>(kTransposed,*this)); }
0150    TMatrixTSparse<Element> &Transpose (const TMatrixTSparse<Element> &source);
0151    inline TMatrixTSparse<Element> &T () { return this->Transpose(*this); }
0152 
0153    inline void Mult(const TMatrixTSparse<Element> &a,const TMatrixTSparse<Element> &b) { AMultB(a,b,0); }
0154 
0155    TMatrixTBase<Element> &Zero       () override;
0156    TMatrixTBase<Element> &UnitMatrix () override;
0157 
0158    Element RowNorm () const override;
0159    Element ColNorm () const override;
0160    Int_t   NonZeros() const override { return this->fNelems; }
0161 
0162    TMatrixTBase<Element> &NormByDiag(const TVectorT<Element> &/*v*/,Option_t * /*option*/) override
0163                                               { MayNotUse("NormByDiag"); return *this; }
0164 
0165    // Either access a_ij as a(i,j)
0166    Element  operator()(Int_t rown,Int_t coln) const override;
0167    Element &operator()(Int_t rown,Int_t coln) override;
0168 
0169    // or as a[i][j]
0170    inline const TMatrixTSparseRow_const<Element> operator[](Int_t rown) const { return TMatrixTSparseRow_const<Element>(*this,rown); }
0171    inline       TMatrixTSparseRow      <Element> operator[](Int_t rown)       { return TMatrixTSparseRow      <Element>(*this,rown); }
0172 
0173    TMatrixTSparse<Element> &operator=(const TMatrixT<Element>       &source);
0174    TMatrixTSparse<Element> &operator=(const TMatrixTSparse<Element> &source);
0175 
0176    TMatrixTSparse<Element> &operator= (Element val);
0177    TMatrixTSparse<Element> &operator-=(Element val);
0178    TMatrixTSparse<Element> &operator+=(Element val);
0179    TMatrixTSparse<Element> &operator*=(Element val);
0180 
0181    TMatrixTSparse<Element> &operator+=(const TMatrixTSparse<Element> &source) { TMatrixTSparse<Element> tmp(*this); Clear();
0182                                                                                 if (this == &source) APlusB (tmp,tmp,1);
0183                                                                                 else                 APlusB (tmp,source,1);
0184                                                                                 return *this; }
0185    TMatrixTSparse<Element> &operator+=(const TMatrixT<Element>       &source) { TMatrixTSparse<Element> tmp(*this); Clear();
0186                                                                                 APlusB(tmp,source,1); return *this; }
0187    TMatrixTSparse<Element> &operator-=(const TMatrixTSparse<Element> &source) { TMatrixTSparse<Element> tmp(*this); Clear();
0188                                                                                 if (this == &source) AMinusB (tmp,tmp,1);
0189                                                                                 else                 AMinusB(tmp,source,1);
0190                                                                                 return *this; }
0191    TMatrixTSparse<Element> &operator-=(const TMatrixT<Element>       &source) { TMatrixTSparse<Element> tmp(*this); Clear();
0192                                                                                 AMinusB(tmp,source,1); return *this; }
0193    TMatrixTSparse<Element> &operator*=(const TMatrixTSparse<Element> &source) { TMatrixTSparse<Element> tmp(*this); Clear();
0194                                                                                 if (this == &source) AMultB (tmp,tmp,1);
0195                                                                                 else                 AMultB (tmp,source,1);
0196                                                                                 return *this; }
0197    TMatrixTSparse<Element> &operator*=(const TMatrixT<Element>       &source) { TMatrixTSparse<Element> tmp(*this); Clear();
0198                                                                                 AMultB(tmp,source,1);
0199                                                                                 return *this; }
0200 
0201    TMatrixTBase  <Element> &Randomize  (Element alpha,Element beta,Double_t &seed) override;
0202    virtual TMatrixTSparse<Element> &RandomizePD(Element alpha,Element beta,Double_t &seed);
0203 
0204    ClassDefOverride(TMatrixTSparse,3) // Template of Sparse Matrix class
0205 };
0206 
0207 // When building with -fmodules, it instantiates all pending instantiations,
0208 // instead of delaying them until the end of the translation unit.
0209 // We 'got away with' probably because the use and the definition of the
0210 // explicit specialization do not occur in the same TU.
0211 //
0212 // In case we are building with -fmodules, we need to forward declare the
0213 // specialization in order to compile the dictionary G__Matrix.cxx.
0214 template <> TClass *TMatrixTSparse<double>::Class();
0215 
0216 template <class Element> inline const Element *TMatrixTSparse<Element>::GetMatrixArray  () const { return fElements; }
0217 template <class Element> inline       Element *TMatrixTSparse<Element>::GetMatrixArray  ()       { return fElements; }
0218 template <class Element> inline const Int_t   *TMatrixTSparse<Element>::GetRowIndexArray() const { return fRowIndex; }
0219 template <class Element> inline       Int_t   *TMatrixTSparse<Element>::GetRowIndexArray()       { return fRowIndex; }
0220 template <class Element> inline const Int_t   *TMatrixTSparse<Element>::GetColIndexArray() const { return fColIndex; }
0221 template <class Element> inline       Int_t   *TMatrixTSparse<Element>::GetColIndexArray()       { return fColIndex; }
0222 
0223 template <class Element>
0224 inline       TMatrixTSparse<Element> &TMatrixTSparse<Element>::Use   (Int_t nrows,Int_t ncols,Int_t nr_nonzeros,
0225                                                                       Int_t *pRowIndex,Int_t *pColIndex,Element *pData)
0226                                                                         { return Use(0,nrows-1,0,ncols-1,nr_nonzeros,pRowIndex,pColIndex,pData); }
0227 template <class Element>
0228 inline const TMatrixTSparse<Element> &TMatrixTSparse<Element>::Use   (Int_t nrows,Int_t ncols,Int_t nr_nonzeros,
0229                                                                       const Int_t *pRowIndex,const Int_t *pColIndex,const Element *pData) const
0230                                                                         { return Use(0,nrows-1,0,ncols-1,nr_nonzeros,pRowIndex,pColIndex,pData); }
0231 template <class Element>
0232 inline       TMatrixTSparse<Element> &TMatrixTSparse<Element>::Use   (TMatrixTSparse<Element> &a)
0233                                                                         { R__ASSERT(a.IsValid());
0234                                                                            return Use(a.GetRowLwb(),a.GetRowUpb(),a.GetColLwb(),a.GetColUpb(),
0235                                                                                       a.GetNoElements(),a.GetRowIndexArray(),
0236                                                                                       a.GetColIndexArray(),a.GetMatrixArray()); }
0237 template <class Element>
0238 inline const TMatrixTSparse<Element> &TMatrixTSparse<Element>::Use   (const TMatrixTSparse<Element> &a) const
0239                                                                         { R__ASSERT(a.IsValid());
0240                                                                            return Use(a.GetRowLwb(),a.GetRowUpb(),a.GetColLwb(),a.GetColUpb(),
0241                                                                                       a.GetNoElements(),a.GetRowIndexArray(),
0242                                                                                       a.GetColIndexArray(),a.GetMatrixArray()); }
0243 
0244 template <class Element>
0245 inline       TMatrixTSparse<Element>  TMatrixTSparse<Element>::GetSub(Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb,
0246                                                                       Option_t *option) const
0247                                                                         {
0248                                                                           TMatrixTSparse<Element> tmp;
0249                                                                           this->GetSub(row_lwb,row_upb,col_lwb,col_upb,tmp,option);
0250                                                                           return tmp;
0251                                                                         }
0252 
0253 template <class Element> TMatrixTSparse<Element>  operator+ (const TMatrixTSparse<Element> &source1,const TMatrixTSparse<Element> &source2);
0254 template <class Element> TMatrixTSparse<Element>  operator+ (const TMatrixTSparse<Element> &source1,const TMatrixT<Element>       &source2);
0255 template <class Element> TMatrixTSparse<Element>  operator+ (const TMatrixT<Element>       &source1,const TMatrixTSparse<Element> &source2);
0256 template <class Element> TMatrixTSparse<Element>  operator+ (const TMatrixTSparse<Element> &source ,      Element                  val    );
0257 template <class Element> TMatrixTSparse<Element>  operator+ (      Element                  val    ,const TMatrixTSparse<Element> &source );
0258 template <class Element> TMatrixTSparse<Element>  operator- (const TMatrixTSparse<Element> &source1,const TMatrixTSparse<Element> &source2);
0259 template <class Element> TMatrixTSparse<Element>  operator- (const TMatrixTSparse<Element> &source1,const TMatrixT<Element>       &source2);
0260 template <class Element> TMatrixTSparse<Element>  operator- (const TMatrixT<Element>       &source1,const TMatrixTSparse<Element> &source2);
0261 template <class Element> TMatrixTSparse<Element>  operator- (const TMatrixTSparse<Element> &source ,      Element                  val    );
0262 template <class Element> TMatrixTSparse<Element>  operator- (      Element                  val    ,const TMatrixTSparse<Element> &source );
0263 template <class Element> TMatrixTSparse<Element>  operator* (const TMatrixTSparse<Element> &source1,const TMatrixTSparse<Element> &source2);
0264 template <class Element> TMatrixTSparse<Element>  operator* (const TMatrixTSparse<Element> &source1,const TMatrixT<Element>       &source2);
0265 template <class Element> TMatrixTSparse<Element>  operator* (const TMatrixT<Element>       &source1,const TMatrixTSparse<Element> &source2);
0266 template <class Element> TMatrixTSparse<Element>  operator* (      Element                  val    ,const TMatrixTSparse<Element> &source );
0267 template <class Element> TMatrixTSparse<Element>  operator* (const TMatrixTSparse<Element> &source,       Element                  val    );
0268 
0269 template <class Element> TMatrixTSparse<Element> &Add        (TMatrixTSparse<Element> &target,      Element                   scalar,
0270                                                               const TMatrixTSparse<Element> &source);
0271 template <class Element> TMatrixTSparse<Element> &ElementMult(TMatrixTSparse<Element> &target,const TMatrixTSparse<Element>  &source);
0272 template <class Element> TMatrixTSparse<Element> &ElementDiv (TMatrixTSparse<Element> &target,const TMatrixTSparse<Element>  &source);
0273 
0274 template <class Element> Bool_t AreCompatible(const TMatrixTSparse<Element> &m1,const TMatrixTSparse<Element> &m2,Int_t verbose=0);
0275 
0276 #endif
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!