Logo ROOT   6.08/07
Reference Guide
TMatrixTUtils.h
Go to the documentation of this file.
1 // @(#)root/matrix:$Id$
2 // Authors: Fons Rademakers, Eddy Offermann Nov 2003
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 #ifndef ROOT_TMatrixTUtils
13 #define ROOT_TMatrixTUtils
14 
15 //////////////////////////////////////////////////////////////////////////
16 // //
17 // Matrix utility classes. //
18 // //
19 // Templates of utility classes in the Linear Algebra Package. //
20 // The following classes are defined here: //
21 // //
22 // Different matrix views without copying data elements : //
23 // TMatrixTRow_const TMatrixTRow //
24 // TMatrixTColumn_const TMatrixTColumn //
25 // TMatrixTDiag_const TMatrixTDiag //
26 // TMatrixTFlat_const TMatrixTFlat //
27 // TMatrixTSub_const TMatrixTSub //
28 // TMatrixTSparseRow_const TMatrixTSparseRow //
29 // TMatrixTSparseDiag_const TMatrixTSparseDiag //
30 // //
31 // TElementActionT //
32 // TElementPosActionT //
33 // //
34 //////////////////////////////////////////////////////////////////////////
35 
36 #ifndef ROOT_TMatrixTBase
37 #include "TMatrixTBase.h"
38 #endif
39 
40 template<class Element> class TVectorT;
41 template<class Element> class TMatrixT;
42 template<class Element> class TMatrixTSym;
43 template<class Element> class TMatrixTSparse;
44 
45 //////////////////////////////////////////////////////////////////////////
46 // //
47 // TElementActionT //
48 // //
49 // A class to do a specific operation on every vector or matrix element //
50 // (regardless of it position) as the object is being traversed. //
51 // This is an abstract class. Derived classes need to implement the //
52 // action function Operation(). //
53 // //
54 //////////////////////////////////////////////////////////////////////////
55 
56 template<class Element> class TElementActionT {
57 
58 #ifndef __CINT__
59 friend class TMatrixTBase <Element>;
60 friend class TMatrixT <Element>;
61 friend class TMatrixTSym <Element>;
62 friend class TMatrixTSparse<Element>;
63 friend class TVectorT <Element>;
64 #endif
65 
66 protected:
67  virtual ~TElementActionT() { }
68  virtual void Operation(Element &element) const = 0;
69 
70 private:
72 };
73 
74 //////////////////////////////////////////////////////////////////////////
75 // //
76 // TElementPosActionT //
77 // //
78 // A class to do a specific operation on every vector or matrix element //
79 // as the object is being traversed. This is an abstract class. //
80 // Derived classes need to implement the action function Operation(). //
81 // In the action function the location of the current element is //
82 // known (fI=row, fJ=columns). //
83 // //
84 //////////////////////////////////////////////////////////////////////////
85 
86 template<class Element> class TElementPosActionT {
87 
88 #ifndef __CINT__
89 friend class TMatrixTBase <Element>;
90 friend class TMatrixT <Element>;
91 friend class TMatrixTSym <Element>;
92 friend class TMatrixTSparse<Element>;
93 friend class TVectorT <Element>;
94 #endif
95 
96 protected:
97  mutable Int_t fI; // i position of element being passed to Operation()
98  mutable Int_t fJ; // j position of element being passed to Operation()
99  virtual ~TElementPosActionT() { }
100  virtual void Operation(Element &element) const = 0;
101 
102 private:
104 };
105 
106 //////////////////////////////////////////////////////////////////////////
107 // //
108 // TMatrixTRow_const //
109 // //
110 // Template class represents a row of a TMatrixT/TMatrixTSym //
111 // //
112 //////////////////////////////////////////////////////////////////////////
113 
114 template<class Element> class TMatrixTRow_const {
115 
116 protected:
117  const TMatrixTBase<Element> *fMatrix; // the matrix I am a row of
118  Int_t fRowInd; // effective row index
119  Int_t fInc; // if ptr = @a[row,i], then ptr+inc = @a[row,i+1]
120  const Element *fPtr; // pointer to the a[row,0]
121 
122 public:
123  TMatrixTRow_const() { fMatrix = 0; fRowInd = 0; fInc = 0; fPtr = 0; }
124  TMatrixTRow_const(const TMatrixT <Element> &matrix,Int_t row);
125  TMatrixTRow_const(const TMatrixTSym<Element> &matrix,Int_t row);
127  fMatrix(trc.fMatrix), fRowInd(trc.fRowInd), fInc(trc.fInc), fPtr(trc.fPtr) { }
129  if(this != &trc) { fMatrix=trc.fMatrix; fRowInd=trc.fRowInd; fInc=trc.fInc; fPtr=trc.fPtr; } return *this;}
130  virtual ~TMatrixTRow_const() { }
131 
132  inline const TMatrixTBase<Element> *GetMatrix () const { return fMatrix; }
133  inline Int_t GetRowIndex() const { return fRowInd; }
134  inline Int_t GetInc () const { return fInc; }
135  inline const Element *GetPtr () const { return fPtr; }
136  inline const Element &operator ()(Int_t i) const {
137  if (!fMatrix) return TMatrixTBase<Element>::NaNValue();
138  R__ASSERT(fMatrix->IsValid());
139  const Int_t acoln = i-fMatrix->GetColLwb();
140  if (acoln < fMatrix->GetNcols() && acoln >= 0)
141  return fPtr[acoln];
142  else {
143  Error("operator()","Request col(%d) outside matrix range of %d - %d",
144  i,fMatrix->GetColLwb(),fMatrix->GetColLwb()+fMatrix->GetNcols());
146  }
147  }
148  inline const Element &operator [](Int_t i) const { return (*(const TMatrixTRow_const<Element> *)this)(i); }
149 
150  ClassDef(TMatrixTRow_const,0) // Template of General Matrix Row Access class
151 };
152 
153 template<class Element> class TMatrixTRow : public TMatrixTRow_const<Element> {
154 
155 public:
157  TMatrixTRow(TMatrixT <Element> &matrix,Int_t row);
160 
161  inline Element *GetPtr() const { return const_cast<Element *>(this->fPtr); }
162 
163  inline const Element &operator()(Int_t i) const {
164  if (!this->fMatrix) return TMatrixTBase<Element>::NaNValue();
165  R__ASSERT(this->fMatrix->IsValid());
166  const Int_t acoln = i-this->fMatrix->GetColLwb();
167  if (acoln < this->fMatrix->GetNcols() || acoln >= 0)
168  return (this->fPtr)[acoln];
169  else {
170  Error("operator()","Request col(%d) outside matrix range of %d - %d",
171  i,this->fMatrix->GetColLwb(),this->fMatrix->GetColLwb()+this->fMatrix->GetNcols());
173  }
174  }
175  inline Element &operator()(Int_t i) {
176  if (!this->fMatrix) return TMatrixTBase<Element>::NaNValue();
177  R__ASSERT(this->fMatrix->IsValid());
178  const Int_t acoln = i-this->fMatrix->GetColLwb();
179  if (acoln < this->fMatrix->GetNcols() && acoln >= 0)
180  return (const_cast<Element *>(this->fPtr))[acoln];
181  else {
182  Error("operator()","Request col(%d) outside matrix range of %d - %d",
183  i,this->fMatrix->GetColLwb(),this->fMatrix->GetColLwb()+this->fMatrix->GetNcols());
184  //return (const_cast<Element *>(this->fPtr))[0];
186  }
187  }
188  inline const Element &operator[](Int_t i) const { return (*(const TMatrixTRow<Element> *)this)(i); }
189  inline Element &operator[](Int_t i) { return (*( TMatrixTRow<Element> *)this)(i); }
190 
191  void operator= (Element val);
192  void operator+=(Element val);
193  void operator*=(Element val);
194 
197  void operator=(const TVectorT <Element> &vec);
198 
199  void operator+=(const TMatrixTRow_const<Element> &r);
200  void operator*=(const TMatrixTRow_const<Element> &r);
201 
202  ClassDef(TMatrixTRow,0) // Template of General Matrix Row Access class
203 };
204 
205 //////////////////////////////////////////////////////////////////////////
206 // //
207 // TMatrixTColumn_const //
208 // //
209 // Template class represents a column of a TMatrixT/TMatrixTSym //
210 // //
211 //////////////////////////////////////////////////////////////////////////
212 
213 template<class Element> class TMatrixTColumn_const {
214 
215 protected:
216  const TMatrixTBase<Element> *fMatrix; // the matrix I am a column of
217  Int_t fColInd; // effective column index
218  Int_t fInc; // if ptr = @a[i,col], then ptr+inc = @a[i+1,col]
219  const Element *fPtr; // pointer to the a[0,col] column
220 
221 public:
222  TMatrixTColumn_const() { fMatrix = 0; fColInd = 0; fInc = 0; fPtr = 0; }
223  TMatrixTColumn_const(const TMatrixT <Element> &matrix,Int_t col);
226  fMatrix(trc.fMatrix), fColInd(trc.fColInd), fInc(trc.fInc), fPtr(trc.fPtr) { }
228  if(this != &trc) { fMatrix=trc.fMatrix; fColInd=trc.fColInd; fInc=trc.fInc; fPtr=trc.fPtr; } return *this;}
229  virtual ~TMatrixTColumn_const() { }
230 
231  inline const TMatrixTBase <Element> *GetMatrix () const { return fMatrix; }
232  inline Int_t GetColIndex() const { return fColInd; }
233  inline Int_t GetInc () const { return fInc; }
234  inline const Element *GetPtr () const { return fPtr; }
235  inline const Element &operator ()(Int_t i) const {
236  if (!this->fMatrix) return TMatrixTBase<Element>::NaNValue();
237  R__ASSERT(fMatrix->IsValid());
238  const Int_t arown = i-fMatrix->GetRowLwb();
239  if (arown < fMatrix->GetNrows() && arown >= 0)
240  return fPtr[arown*fInc];
241  else {
242  Error("operator()","Request row(%d) outside matrix range of %d - %d",
243  i,fMatrix->GetRowLwb(),fMatrix->GetRowLwb()+fMatrix->GetNrows());
245  }
246  }
247  inline const Element &operator [](Int_t i) const { return (*(const TMatrixTColumn_const<Element> *)this)(i); }
248 
249  ClassDef(TMatrixTColumn_const,0) // Template of General Matrix Column Access class
250 };
251 
252 template<class Element> class TMatrixTColumn : public TMatrixTColumn_const<Element> {
253 
254 public:
259 
260  inline Element *GetPtr() const { return const_cast<Element *>(this->fPtr); }
261 
262  inline const Element &operator()(Int_t i) const {
263  if (!this->fMatrix) return TMatrixTBase<Element>::NaNValue();
264  R__ASSERT(this->fMatrix->IsValid());
265  const Int_t arown = i-this->fMatrix->GetRowLwb();
266  if (arown < this->fMatrix->GetNrows() && arown >= 0)
267  return (this->fPtr)[arown*this->fInc];
268  else {
269  Error("operator()","Request row(%d) outside matrix range of %d - %d",
270  i,this->fMatrix->GetRowLwb(),this->fMatrix->GetRowLwb()+this->fMatrix->GetNrows());
272  }
273  }
274  inline Element &operator()(Int_t i) {
275  if (!this->fMatrix) return TMatrixTBase<Element>::NaNValue();
276  R__ASSERT(this->fMatrix->IsValid());
277  const Int_t arown = i-this->fMatrix->GetRowLwb();
278 
279  if (arown < this->fMatrix->GetNrows() && arown >= 0)
280  return (const_cast<Element *>(this->fPtr))[arown*this->fInc];
281  else {
282  Error("operator()","Request row(%d) outside matrix range of %d - %d",
283  i,this->fMatrix->GetRowLwb(),this->fMatrix->GetRowLwb()+this->fMatrix->GetNrows());
285  }
286  }
287  inline const Element &operator[](Int_t i) const { return (*(const TMatrixTColumn<Element> *)this)(i); }
288  inline Element &operator[](Int_t i) { return (*( TMatrixTColumn<Element> *)this)(i); }
289 
290  void operator= (Element val);
291  void operator+=(Element val);
292  void operator*=(Element val);
293 
296  void operator=(const TVectorT <Element> &vec);
297 
298  void operator+=(const TMatrixTColumn_const<Element> &c);
299  void operator*=(const TMatrixTColumn_const<Element> &c);
300 
301  ClassDef(TMatrixTColumn,0) // Template of General Matrix Column Access class
302 };
303 
304 //////////////////////////////////////////////////////////////////////////
305 // //
306 // TMatrixTDiag_const //
307 // //
308 // Template class represents the diagonal of a TMatrixT/TMatrixTSym //
309 // //
310 //////////////////////////////////////////////////////////////////////////
311 
312 template<class Element> class TMatrixTDiag_const {
313 
314 protected:
315  const TMatrixTBase<Element> *fMatrix; // the matrix I am the diagonal of
316  Int_t fInc; // if ptr=@a[i,i], then ptr+inc = @a[i+1,i+1]
317  Int_t fNdiag; // number of diag elems, min(nrows,ncols)
318  const Element *fPtr; // pointer to the a[0,0]
319 
320 public:
321  TMatrixTDiag_const() { fMatrix = 0; fInc = 0; fNdiag = 0; fPtr = 0; }
322  TMatrixTDiag_const(const TMatrixT <Element> &matrix);
325  fMatrix(trc.fMatrix), fInc(trc.fInc), fNdiag(trc.fNdiag), fPtr(trc.fPtr) { }
327  if(this != &trc) { fMatrix=trc.fMatrix; fInc=trc.fInc; fNdiag=trc.fNdiag; fPtr=trc.fPtr; } return *this;}
328  virtual ~TMatrixTDiag_const() { }
329 
330  inline const TMatrixTBase<Element> *GetMatrix() const { return fMatrix; }
331  inline const Element *GetPtr () const { return fPtr; }
332  inline Int_t GetInc () const { return fInc; }
333  inline const Element &operator ()(Int_t i) const {
334  R__ASSERT(fMatrix->IsValid());
335  if (i < fNdiag && i >= 0)
336  return fPtr[i*fInc];
337  else {
338  Error("operator()","Request diagonal(%d) outside matrix range of 0 - %d",i,fNdiag);
340  }
341  }
342  inline const Element &operator [](Int_t i) const { return (*(const TMatrixTDiag_const<Element> *)this)(i); }
343 
344  Int_t GetNdiags() const { return fNdiag; }
345 
346  ClassDef(TMatrixTDiag_const,0) // Template of General Matrix Diagonal Access class
347 };
348 
349 template<class Element> class TMatrixTDiag : public TMatrixTDiag_const<Element> {
350 
351 public:
356 
357  inline Element *GetPtr() const { return const_cast<Element *>(this->fPtr); }
358 
359  inline const Element &operator()(Int_t i) const {
360  R__ASSERT(this->fMatrix->IsValid());
361  if (i < this->fNdiag && i >= 0)
362  return (this->fPtr)[i*this->fInc];
363  else {
364  Error("operator()","Request diagonal(%d) outside matrix range of 0 - %d",i,this->fNdiag);
366  }
367  }
368  inline Element &operator()(Int_t i) {
369  R__ASSERT(this->fMatrix->IsValid());
370  if (i < this->fNdiag && i >= 0)
371  return (const_cast<Element *>(this->fPtr))[i*this->fInc];
372  else {
373  Error("operator()","Request diagonal(%d) outside matrix range of 0 - %d",i,this->fNdiag);
374  return (const_cast<Element *>(this->fPtr))[0];
375  }
376  }
377  inline const Element &operator[](Int_t i) const { return (*(const TMatrixTDiag<Element> *)this)(i); }
378  inline Element &operator[](Int_t i) { return (*( TMatrixTDiag *)this)(i); }
379 
380  void operator= (Element val);
381  void operator+=(Element val);
382  void operator*=(Element val);
383 
384  void operator=(const TMatrixTDiag_const<Element> &d);
386  void operator=(const TVectorT <Element> &vec);
387 
388  void operator+=(const TMatrixTDiag_const<Element> &d);
389  void operator*=(const TMatrixTDiag_const<Element> &d);
390 
391  ClassDef(TMatrixTDiag,0) // Template of General Matrix Diagonal Access class
392 };
393 
394 //////////////////////////////////////////////////////////////////////////
395 // //
396 // TMatrixTFlat_const //
397 // //
398 // Template class represents a flat TMatrixT/TMatrixTSym //
399 // //
400 //////////////////////////////////////////////////////////////////////////
401 
402 template<class Element> class TMatrixTFlat_const {
403 
404 protected:
405  const TMatrixTBase<Element> *fMatrix; // the matrix I am the diagonal of
407  const Element *fPtr; // pointer to the a[0,0]
408 
409 public:
410  TMatrixTFlat_const() { fMatrix = 0; fNelems = 0; fPtr = 0; }
411  TMatrixTFlat_const(const TMatrixT <Element> &matrix);
414  fMatrix(trc.fMatrix), fNelems(trc.fNelems), fPtr(trc.fPtr) { }
416  if(this != &trc) { fMatrix=trc.fMatrix; fNelems=trc.fNelems; fPtr=trc.fPtr; } return *this;}
417  virtual ~TMatrixTFlat_const() { }
418 
419  inline const TMatrixTBase<Element> *GetMatrix() const { return fMatrix; }
420  inline const Element *GetPtr () const { return fPtr; }
421  inline const Element &operator ()(Int_t i) const {
422  R__ASSERT(fMatrix->IsValid());
423  if (i < fNelems && i >= 0)
424  return fPtr[i];
425  else {
426  Error("operator()","Request element(%d) outside matrix range of 0 - %d",i,fNelems);
428  }
429  }
430  inline const Element &operator [](Int_t i) const { return (*(const TMatrixTFlat_const<Element> *)this)(i); }
431 
432  ClassDef(TMatrixTFlat_const,0) // Template of General Matrix Flat Representation class
433 };
434 
435 template<class Element> class TMatrixTFlat : public TMatrixTFlat_const<Element> {
436 
437 public:
442 
443  inline Element *GetPtr() const { return const_cast<Element *>(this->fPtr); }
444 
445  inline const Element &operator()(Int_t i) const {
446  R__ASSERT(this->fMatrix->IsValid());
447  if (i < this->fNelems && i >= 0)
448  return (this->fPtr)[i];
449  else {
450  Error("operator()","Request element(%d) outside matrix range of 0 - %d",i,this->fNelems);
452  }
453  }
454  inline Element &operator()(Int_t i) {
455  R__ASSERT(this->fMatrix->IsValid());
456  if (i < this->fNelems && i >= 0)
457  return (const_cast<Element *>(this->fPtr))[i];
458  else {
459  Error("operator()","Request element(%d) outside matrix range of 0 - %d",i,this->fNelems);
461  }
462  }
463  inline const Element &operator[](Int_t i) const { return (*(const TMatrixTFlat<Element> *)this)(i); }
464  inline Element &operator[](Int_t i) { return (*( TMatrixTFlat<Element> *)this)(i); }
465 
466  void operator= (Element val);
467  void operator+=(Element val);
468  void operator*=(Element val);
469 
472  void operator=(const TVectorT <Element> &vec);
473 
474  void operator+=(const TMatrixTFlat_const<Element> &f);
475  void operator*=(const TMatrixTFlat_const<Element> &f);
476 
477  ClassDef(TMatrixTFlat,0) // Template of General Matrix Flat Representation class
478 };
479 
480 //////////////////////////////////////////////////////////////////////////
481 // //
482 // TMatrixTSub_const //
483 // //
484 // Template class represents a sub matrix of TMatrixT/TMatrixTSym //
485 // //
486 //////////////////////////////////////////////////////////////////////////
487 
488 template<class Element> class TMatrixTSub_const {
489 
490 protected:
491  const TMatrixTBase<Element> *fMatrix; // the matrix I am a submatrix of
496 
497 public:
498  TMatrixTSub_const() { fRowOff = fColOff = fNrowsSub = fNcolsSub = 0; fMatrix = 0; }
499  TMatrixTSub_const(const TMatrixT <Element> &matrix,Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb);
500  TMatrixTSub_const(const TMatrixTSym<Element> &matrix,Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb);
501  virtual ~TMatrixTSub_const() { }
502 
503  inline const TMatrixTBase<Element> *GetMatrix() const { return fMatrix; }
504  inline Int_t GetRowOff() const { return fRowOff; }
505  inline Int_t GetColOff() const { return fColOff; }
506  inline Int_t GetNrows () const { return fNrowsSub; }
507  inline Int_t GetNcols () const { return fNcolsSub; }
508  inline const Element &operator ()(Int_t rown,Int_t coln) const {
509  R__ASSERT(fMatrix->IsValid());
510 
511  const Element *ptr = fMatrix->GetMatrixArray();
512  if (rown >= fNrowsSub || rown < 0) {
513  Error("operator()","Request row(%d) outside matrix range of 0 - %d",rown,fNrowsSub);
515  }
516  if (coln >= fNcolsSub || coln < 0) {
517  Error("operator()","Request column(%d) outside matrix range of 0 - %d",coln,fNcolsSub);
519  }
520  const Int_t index = (rown+fRowOff)*fMatrix->GetNcols()+coln+fColOff;
521  return ptr[index];
522  }
523 
524  ClassDef(TMatrixTSub_const,0) // Template of Sub Matrix Access class
525 };
526 
527 template<class Element> class TMatrixTSub : public TMatrixTSub_const<Element> {
528 
529 public:
530 
531  enum {kWorkMax = 100};
532 
534  TMatrixTSub(TMatrixT <Element> &matrix,Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb);
535  TMatrixTSub(TMatrixTSym<Element> &matrix,Int_t row_lwb,Int_t row_upb,Int_t col_lwb,Int_t col_upb);
537 
538  inline Element &operator()(Int_t rown,Int_t coln) {
539  R__ASSERT(this->fMatrix->IsValid());
540 
541  const Element *ptr = this->fMatrix->GetMatrixArray();
542  if (rown >= this->fNrowsSub || rown < 0) {
543  Error("operator()","Request row(%d) outside matrix range of 0 - %d",rown,this->fNrowsSub);
545  }
546  if (coln >= this->fNcolsSub || coln < 0) {
547  Error("operator()","Request column(%d) outside matrix range of 0 - %d",coln,this->fNcolsSub);
549  }
550  const Int_t index = (rown+this->fRowOff)*this->fMatrix->GetNcols()+coln+this->fColOff;
551  return (const_cast<Element *>(ptr))[index];
552  }
553 
554  void Rank1Update(const TVectorT<Element> &vec,Element alpha=1.0);
555 
556  void operator= (Element val);
557  void operator+=(Element val);
558  void operator*=(Element val);
559 
560  void operator=(const TMatrixTSub_const<Element> &s);
562  void operator=(const TMatrixTBase <Element> &m);
563 
564  void operator+=(const TMatrixTSub_const<Element> &s);
565  void operator*=(const TMatrixTSub_const<Element> &s);
566  void operator+=(const TMatrixTBase <Element> &m);
567  void operator*=(const TMatrixT <Element> &m);
568  void operator*=(const TMatrixTSym <Element> &m);
569 
570  ClassDef(TMatrixTSub,0) // Template of Sub Matrix Access class
571 };
572 
573 //////////////////////////////////////////////////////////////////////////
574 // //
575 // TMatrixTSparseRow_const //
576 // //
577 // Template class represents a row of TMatrixTSparse //
578 // //
579 //////////////////////////////////////////////////////////////////////////
580 
581 template<class Element> class TMatrixTSparseRow_const {
582 
583 protected:
584  const TMatrixTSparse<Element> *fMatrix; // the matrix I am a row of
585  Int_t fRowInd; // effective row index
586  Int_t fNindex; // index range
587  const Int_t *fColPtr; // column index pointer
588  const Element *fDataPtr; // data pointer
589 
590 public:
591  TMatrixTSparseRow_const() { fMatrix = 0; fRowInd = 0; fNindex = 0; fColPtr = 0; fDataPtr = 0; }
594  fMatrix(trc.fMatrix), fRowInd(trc.fRowInd), fNindex(trc.fNindex), fColPtr(trc.fColPtr), fDataPtr(trc.fDataPtr) { }
596  if(this != &trc) { fMatrix=trc.fMatrix; fRowInd=trc.fRowInd; fNindex=trc.fNindex; fColPtr=trc.fColPtr; fDataPtr=trc.fDataPtr; } return *this;}
598 
599  inline const TMatrixTBase<Element> *GetMatrix () const { return fMatrix; }
600  inline const Element *GetDataPtr () const { return fDataPtr; }
601  inline const Int_t *GetColPtr () const { return fColPtr; }
602  inline Int_t GetRowIndex() const { return fRowInd; }
603  inline Int_t GetNindex () const { return fNindex; }
604 
605  Element operator()(Int_t i) const;
606  inline Element operator[](Int_t i) const { return (*(const TMatrixTSparseRow_const<Element> *)this)(i); }
607 
608  ClassDef(TMatrixTSparseRow_const,0) // Template of Sparse Matrix Row Access class
609 };
610 
611 template<class Element> class TMatrixTSparseRow : public TMatrixTSparseRow_const<Element> {
612 
613 public:
617 
618  inline Element *GetDataPtr() const { return const_cast<Element *>(this->fDataPtr); }
619 
620  Element operator()(Int_t i) const;
621  Element &operator()(Int_t i);
622  inline Element operator[](Int_t i) const { return (*(const TMatrixTSparseRow<Element> *)this)(i); }
623  inline Element &operator[](Int_t i) { return (*(TMatrixTSparseRow<Element> *)this)(i); }
624 
625  void operator= (Element val);
626  void operator+=(Element val);
627  void operator*=(Element val);
628 
631  void operator=(const TVectorT <Element> &vec);
632 
633  void operator+=(const TMatrixTSparseRow_const<Element> &r);
634  void operator*=(const TMatrixTSparseRow_const<Element> &r);
635 
636  ClassDef(TMatrixTSparseRow,0) // Template of Sparse Matrix Row Access class
637 };
638 
639 //////////////////////////////////////////////////////////////////////////
640 // //
641 // TMatrixTSparseDiag_const //
642 // //
643 // Template class represents the diagonal of TMatrixTSparse //
644 // //
645 //////////////////////////////////////////////////////////////////////////
646 
647 template<class Element> class TMatrixTSparseDiag_const {
648 
649 protected:
650  const TMatrixTSparse<Element> *fMatrix; // the matrix I am the diagonal of
651  Int_t fNdiag; // number of diag elems, min(nrows,ncols)
652  const Element *fDataPtr; // data pointer
653 
654 public:
655  TMatrixTSparseDiag_const() { fMatrix = 0; fNdiag = 0; fDataPtr = 0; }
658  fMatrix(trc.fMatrix), fNdiag(trc.fNdiag), fDataPtr(trc.fDataPtr) { }
660  if(this != &trc) { fMatrix=trc.fMatrix; fNdiag=trc.fNdiag; fDataPtr=trc.fDataPtr; } return *this;}
662 
663  inline const TMatrixTBase<Element> *GetMatrix () const { return fMatrix; }
664  inline const Element *GetDataPtr() const { return fDataPtr; }
665  inline Int_t GetNdiags () const { return fNdiag; }
666 
667  Element operator ()(Int_t i) const;
668  inline Element operator [](Int_t i) const { return (*(const TMatrixTSparseRow_const<Element> *)this)(i); }
669 
670  ClassDef(TMatrixTSparseDiag_const,0) // Template of Sparse Matrix Diagonal Access class
671 };
672 
673 template<class Element> class TMatrixTSparseDiag : public TMatrixTSparseDiag_const<Element> {
674 
675 public:
679 
680  inline Element *GetDataPtr() const { return const_cast<Element *>(this->fDataPtr); }
681 
682  Element operator()(Int_t i) const;
683  Element &operator()(Int_t i);
684  inline Element operator[](Int_t i) const { return (*(const TMatrixTSparseDiag<Element> *)this)(i); }
685  inline Element &operator[](Int_t i) { return (*(TMatrixTSparseDiag<Element> *)this)(i); }
686 
687  void operator= (Element val);
688  void operator+=(Element val);
689  void operator*=(Element val);
690 
693  void operator=(const TVectorT <Element> &vec);
694 
695  void operator+=(const TMatrixTSparseDiag_const<Element> &d);
696  void operator*=(const TMatrixTSparseDiag_const<Element> &d);
697 
698  ClassDef(TMatrixTSparseDiag,0) // Template of Sparse Matrix Diagonal Access class
699 };
700 
701 Double_t Drand(Double_t &ix);
702 #endif
const TMatrixTSparse< Element > * fMatrix
TMatrixTSparseDiag_const< Element > & operator=(const TMatrixTSparseDiag_const< Element > &trc)
virtual const Element * GetMatrixArray() const =0
const TMatrixTBase< Element > * GetMatrix() const
const Element * GetPtr() const
Element & operator()(Int_t i)
Int_t GetRowLwb() const
Definition: TMatrixTBase.h:132
Int_t GetInc() const
Int_t GetRowIndex() const
return c
const Element * fPtr
virtual ~TMatrixTSparseDiag_const()
TMatrixTSparseRow< Element > & operator=(const TMatrixTSparseRow< Element > &r)
const TMatrixTBase< Element > * fMatrix
TMatrixTFlat_const(const TMatrixTFlat_const< Element > &trc)
Int_t GetNcols() const
Definition: TMatrixTBase.h:137
TVectorT.
Definition: TMatrixTBase.h:89
#define R__ASSERT(e)
Definition: TError.h:98
Element & operator()(Int_t i)
TMatrixTColumn< Element > & operator=(const TMatrixTColumn< Element > &c)
int Int_t
Definition: RtypesCore.h:41
const Element * fDataPtr
const Element & operator()(Int_t i) const
Element & operator()(Int_t i)
const Element & operator[](Int_t i) const
Int_t GetInc() const
const Element * GetPtr() const
TMatrixTSparseRow_const< Element > & operator=(const TMatrixTSparseRow_const< Element > &trc)
TMatrixTSub< Element > & operator=(const TMatrixTSub< Element > &s)
Element operator[](Int_t i) const
Element & operator()(Int_t rown, Int_t coln)
TMatrixT.
Definition: TMatrixDfwd.h:24
TMatrixTDiag_const< Element > & operator=(const TMatrixTDiag_const< Element > &trc)
virtual ~TElementPosActionT()
Definition: TMatrixTUtils.h:99
Int_t GetNrows() const
const Element * GetPtr() const
const Element * fPtr
const Element & operator()(Int_t i) const
#define ClassDef(name, id)
Definition: Rtypes.h:254
Element * GetPtr() const
Int_t GetColLwb() const
Definition: TMatrixTBase.h:135
Element * GetPtr() const
virtual void Operation(Element &element) const =0
TMatrixTSparseDiag_const(const TMatrixTSparseDiag_const< Element > &trc)
const TMatrixTBase< Element > * fMatrix
Element & operator[](Int_t i)
const Element * GetDataPtr() const
const Element * GetPtr() const
Element * GetDataPtr() const
Element operator[](Int_t i) const
TMatrixTSparseDiag< Element > & operator=(const TMatrixTSparseDiag< Element > &d)
const Element * fDataPtr
Int_t GetRowIndex() const
TMatrixTSym.
const Element * GetDataPtr() const
const Element & operator[](Int_t i) const
void Error(const char *location, const char *msgfmt,...)
Int_t GetColOff() const
virtual ~TMatrixTSub_const()
Element * GetDataPtr() const
TMatrixTDiag< Element > & operator=(const TMatrixTDiag< Element > &d)
Element & operator[](Int_t i)
const Element & operator()(Int_t i) const
TRandom2 r(17)
TMatrixTColumn_const< Element > & operator=(const TMatrixTColumn_const< Element > &trc)
TMatrixTSparseRow_const(const TMatrixTSparseRow_const< Element > &trc)
TMatrixTSparse.
Int_t GetNcols() const
const Element & operator()(Int_t i) const
TMarker * m
Definition: textangle.C:8
const Element * fPtr
Element & operator[](Int_t i)
const TMatrixTBase< Element > * fMatrix
TMatrixTRow< Element > & operator=(const TMatrixTRow< Element > &r)
TElementActionT & operator=(const TElementActionT< Element > &)
Definition: TMatrixTUtils.h:71
Linear Algebra Package.
const TMatrixTBase< Element > * GetMatrix() const
virtual ~TMatrixTRow_const()
Element operator[](Int_t i) const
const Int_t * GetColPtr() const
Int_t GetNrows() const
Definition: TMatrixTBase.h:134
virtual ~TMatrixTFlat_const()
Element * GetPtr() const
const TMatrixTBase< Element > * GetMatrix() const
Element * GetPtr() const
Element & operator[](Int_t i)
TMatrixTRow_const< Element > & operator=(const TMatrixTRow_const< Element > &trc)
double f(double x)
double Double_t
Definition: RtypesCore.h:55
TMatrixTRow_const(const TMatrixTRow_const< Element > &trc)
const TMatrixTBase< Element > * GetMatrix() const
const TMatrixTSparse< Element > * fMatrix
virtual ~TElementActionT()
Definition: TMatrixTUtils.h:67
const TMatrixTBase< Element > * fMatrix
Int_t GetRowOff() const
const TMatrixTBase< Element > * GetMatrix() const
TMatrixTDiag_const(const TMatrixTDiag_const< Element > &trc)
Element & operator[](Int_t i)
const Element * fPtr
Int_t GetColIndex() const
const Element & operator[](Int_t i) const
virtual ~TMatrixTDiag_const()
static Element & NaNValue()
Element & operator()(Int_t i)
virtual ~TMatrixTColumn_const()
const TMatrixTBase< Element > * fMatrix
TMatrixTColumn_const(const TMatrixTColumn_const< Element > &trc)
virtual ~TMatrixTSparseRow_const()
Double_t Drand(Double_t &ix)
Random number generator [0....1] with seed ix.
const Element & operator[](Int_t i) const
TMatrixTFlat_const< Element > & operator=(const TMatrixTFlat_const< Element > &trc)
const TMatrixTBase< Element > * GetMatrix() const
Element & operator[](Int_t i)
const TMatrixTBase< Element > * GetMatrix() const
Bool_t IsValid() const
Definition: TMatrixTBase.h:157
TMatrixTFlat< Element > & operator=(const TMatrixTFlat< Element > &f)
Int_t GetNdiags() const
Int_t GetInc() const
TElementPosActionT< Element > & operator=(const TElementPosActionT< Element > &)