Logo ROOT  
Reference Guide
MatrixRepresentationsStatic.h
Go to the documentation of this file.
1 // @(#)root/smatrix:$Id$
2 // Author: L. Moneta, J. Palacios 2006
3 
4 #ifndef ROOT_Math_MatrixRepresentationsStatic
5 #define ROOT_Math_MatrixRepresentationsStatic 1
6 
7 // Include files
8 
9 /**
10 \defgroup MatRep SMatrix Storage Representation
11 \ingroup SMatrixGroup
12 
13 Classes MatRepStd and MatRepSym for generic and symmetric matrix
14 data storage and manipulation. Define data storage and access, plus
15 operators =, +=, -=, ==.
16 
17 \author Juan Palacios
18 \date 2006-01-15
19  */
20 
21 #include "Math/StaticCheck.h"
22 
23 #include <cstddef>
24 #include <utility>
25 #include <type_traits>
26 #include <array>
27 
28 namespace ROOT {
29 
30 namespace Math {
31 
32 /**
33 \defgroup MatRepStd Standard Matrix representation
34 \ingroup MatRep
35 
36 Standard Matrix representation for a general D1 x D2 matrix.
37 This class is itself a template on the contained type T, the number of rows and the number of columns.
38 Its data member is an array T[nrows*ncols] containing the matrix data.
39 The data are stored in the row-major C convention.
40 For example, for a matrix, M, of size 3x3, the data \f$ \left[a_0,a_1,a_2,.......,a_7,a_8 \right] \f$d
41 are stored in the following order:
42 
43 \f[
44 M = \left( \begin{array}{ccc}
45 a_0 & a_1 & a_2 \\
46 a_3 & a_4 & a_5 \\
47 a_6 & a_7 & a_8 \end{array} \right)
48 \f]
49 
50 */
51 
52 
53  template <class T, unsigned int D1, unsigned int D2=D1>
54  class MatRepStd {
55 
56  public:
57 
58  typedef T value_type;
59 
60  inline const T& operator()(unsigned int i, unsigned int j) const {
61  return fArray[i*D2+j];
62  }
63  inline T& operator()(unsigned int i, unsigned int j) {
64  return fArray[i*D2+j];
65  }
66  inline T& operator[](unsigned int i) { return fArray[i]; }
67 
68  inline const T& operator[](unsigned int i) const { return fArray[i]; }
69 
70  inline T apply(unsigned int i) const { return fArray[i]; }
71 
72  inline T* Array() { return fArray; }
73 
74  inline const T* Array() const { return fArray; }
75 
76  template <class R>
77  inline MatRepStd<T, D1, D2>& operator+=(const R& rhs) {
78  for(unsigned int i=0; i<kSize; ++i) fArray[i] += rhs[i];
79  return *this;
80  }
81 
82  template <class R>
83  inline MatRepStd<T, D1, D2>& operator-=(const R& rhs) {
84  for(unsigned int i=0; i<kSize; ++i) fArray[i] -= rhs[i];
85  return *this;
86  }
87 
88  template <class R>
89  inline MatRepStd<T, D1, D2>& operator=(const R& rhs) {
90  for(unsigned int i=0; i<kSize; ++i) fArray[i] = rhs[i];
91  return *this;
92  }
93 
94  template <class R>
95  inline bool operator==(const R& rhs) const {
96  bool rc = true;
97  for(unsigned int i=0; i<kSize; ++i) {
98  rc = rc && (fArray[i] == rhs[i]);
99  }
100  return rc;
101  }
102 
103  enum {
104  /// return no. of matrix rows
105  kRows = D1,
106  /// return no. of matrix columns
107  kCols = D2,
108  /// return no of elements: rows*columns
109  kSize = D1*D2
110  };
111 
112  private:
113  //T __attribute__ ((aligned (16))) fArray[kSize];
115  };
116 
117 
118 // template<unigned int D>
119 // struct Creator {
120 // static const RowOffsets<D> & Offsets() {
121 // static RowOffsets<D> off;
122 // return off;
123 // }
124 
125  /**
126  Static structure to keep the conversion from (i,j) to offsets in the storage data for a
127  symmetric matrix
128  */
129 
130  template<unsigned int D>
131  struct RowOffsets {
132  inline RowOffsets() {
133  int v[D];
134  v[0]=0;
135  for (unsigned int i=1; i<D; ++i)
136  v[i]=v[i-1]+i;
137  for (unsigned int i=0; i<D; ++i) {
138  for (unsigned int j=0; j<=i; ++j)
139  fOff[i*D+j] = v[i]+j;
140  for (unsigned int j=i+1; j<D; ++j)
141  fOff[i*D+j] = v[j]+i ;
142  }
143  }
144  inline int operator()(unsigned int i, unsigned int j) const { return fOff[i*D+j]; }
145  inline int apply(unsigned int i) const { return fOff[i]; }
146  int fOff[D*D];
147  };
148 
149  namespace rowOffsetsUtils {
150 
151  ///////////
152  // Some meta template stuff
153  template<int...> struct indices{};
154 
155  template<int I, class IndexTuple, int N>
157 
158  template<int I, int... Indices, int N>
159  struct make_indices_impl<I, indices<Indices...>, N>
160  {
161  typedef typename make_indices_impl<I + 1, indices<Indices..., I>,
163  };
164 
165  template<int N, int... Indices>
166  struct make_indices_impl<N, indices<Indices...>, N> {
167  typedef indices<Indices...> type;
168  };
169 
170  template<int N>
171  struct make_indices : make_indices_impl<0, indices<>, N> {};
172  // end of stuff
173 
174 
175 
176  template<int I0, class F, int... I>
177  constexpr std::array<decltype(std::declval<F>()(std::declval<int>())), sizeof...(I)>
179  {
180  return std::array<decltype(std::declval<F>()(std::declval<int>())),
181  sizeof...(I)>{{ f(I0 + I)... }};
182  }
183 
184  template<int N, int I0 = 0, class F>
185  constexpr std::array<decltype(std::declval<F>()(std::declval<int>())), N>
186  make(F f) {
187  return do_make<I0>(f, typename make_indices<N>::type());
188  }
189 
190  } // namespace rowOffsetsUtils
191 
192 
193 //_________________________________________________________________________________
194  /**
195  MatRepSym
196  Matrix storage representation for a symmetric matrix of dimension NxN
197  This class is a template on the contained type and on the symmetric matrix size, N.
198  It has as data member an array of type T of size N*(N+1)/2,
199  containing the lower diagonal block of the matrix.
200  The order follows the lower diagonal block, still in a row-major convention.
201  For example for a symmetric 3x3 matrix the order of the 6 elements
202  \f$ \left[a_0,a_1.....a_5 \right]\f$ is:
203  \f[
204  M = \left( \begin{array}{ccc}
205  a_0 & a_1 & a_3 \\
206  a_1 & a_2 & a_4 \\
207  a_3 & a_4 & a_5 \end{array} \right)
208  \f]
209 
210  @ingroup MatRep
211  */
212  template <class T, unsigned int D>
213  class MatRepSym {
214 
215  public:
216 
217  /* constexpr */ inline MatRepSym(){}
218 
219  typedef T value_type;
220 
221 
222  inline T & operator()(unsigned int i, unsigned int j)
223  { return fArray[offset(i, j)]; }
224 
225  inline /* constexpr */ T const & operator()(unsigned int i, unsigned int j) const
226  { return fArray[offset(i, j)]; }
227 
228  inline T& operator[](unsigned int i) {
229  return fArray[off(i)];
230  }
231 
232  inline /* constexpr */ T const & operator[](unsigned int i) const {
233  return fArray[off(i)];
234  }
235 
236  inline /* constexpr */ T apply(unsigned int i) const {
237  return fArray[off(i)];
238  }
239 
240  inline T* Array() { return fArray; }
241 
242  inline const T* Array() const { return fArray; }
243 
244  /**
245  assignment : only symmetric to symmetric allowed
246  */
247  template <class R>
248  inline MatRepSym<T, D>& operator=(const R&) {
249  STATIC_CHECK(0==1,
250  Cannot_assign_general_to_symmetric_matrix_representation);
251  return *this;
252  }
253  inline MatRepSym<T, D>& operator=(const MatRepSym& rhs) {
254  for(unsigned int i=0; i<kSize; ++i) fArray[i] = rhs.Array()[i];
255  return *this;
256  }
257 
258  /**
259  self addition : only symmetric to symmetric allowed
260  */
261  template <class R>
262  inline MatRepSym<T, D>& operator+=(const R&) {
263  STATIC_CHECK(0==1,
264  Cannot_add_general_to_symmetric_matrix_representation);
265  return *this;
266  }
267  inline MatRepSym<T, D>& operator+=(const MatRepSym& rhs) {
268  for(unsigned int i=0; i<kSize; ++i) fArray[i] += rhs.Array()[i];
269  return *this;
270  }
271 
272  /**
273  self subtraction : only symmetric to symmetric allowed
274  */
275  template <class R>
276  inline MatRepSym<T, D>& operator-=(const R&) {
277  STATIC_CHECK(0==1,
278  Cannot_substract_general_to_symmetric_matrix_representation);
279  return *this;
280  }
281  inline MatRepSym<T, D>& operator-=(const MatRepSym& rhs) {
282  for(unsigned int i=0; i<kSize; ++i) fArray[i] -= rhs.Array()[i];
283  return *this;
284  }
285  template <class R>
286  inline bool operator==(const R& rhs) const {
287  bool rc = true;
288  for(unsigned int i=0; i<D*D; ++i) {
289  rc = rc && (operator[](i) == rhs[i]);
290  }
291  return rc;
292  }
293 
294  enum {
295  /// return no. of matrix rows
296  kRows = D,
297  /// return no. of matrix columns
298  kCols = D,
299  /// return no of elements: rows*columns
300  kSize = D*(D+1)/2
301  };
302 
303  static constexpr int off0(int i) { return i==0 ? 0 : off0(i-1)+i;}
304  static constexpr int off2(int i, int j) { return j<i ? off0(i)+j : off0(j)+i; }
305  static constexpr int off1(int i) { return off2(i/D, i%D);}
306 
307  static int off(int i) {
308  static constexpr auto v = rowOffsetsUtils::make<D*D>(off1);
309  return v[i];
310  }
311 
312  static inline constexpr unsigned int
313  offset(unsigned int i, unsigned int j)
314  {
315  //if (j > i) std::swap(i, j);
316  return off(i*D+j);
317  // return (i>j) ? (i * (i+1) / 2) + j : (j * (j+1) / 2) + i;
318  }
319 
320  private:
321  //T __attribute__ ((aligned (16))) fArray[kSize];
323  };
324 
325 
326 
327 } // namespace Math
328 } // namespace ROOT
329 
330 
331 #endif // MATH_MATRIXREPRESENTATIONSSTATIC_H
ROOT::Math::MatRepStd::operator+=
MatRepStd< T, D1, D2 > & operator+=(const R &rhs)
Definition: MatrixRepresentationsStatic.h:77
ROOT::Math::rowOffsetsUtils::make
constexpr std::array< decltype(std::declval< F >)(std::declval< int >))), N > make(F f)
Definition: MatrixRepresentationsStatic.h:186
ROOT::Math::MatRepStd::operator[]
const T & operator[](unsigned int i) const
Definition: MatrixRepresentationsStatic.h:68
ROOT::Math::rowOffsetsUtils::make_indices_impl
Definition: MatrixRepresentationsStatic.h:156
ROOT::Math::MatRepSym::MatRepSym
MatRepSym()
Definition: MatrixRepresentationsStatic.h:217
ROOT::Math::MatRepSym::off0
static constexpr int off0(int i)
Definition: MatrixRepresentationsStatic.h:303
ROOT::Math::MatRepStd::operator()
T & operator()(unsigned int i, unsigned int j)
Definition: MatrixRepresentationsStatic.h:63
ROOT::Math::MatRepSym::operator+=
MatRepSym< T, D > & operator+=(const MatRepSym &rhs)
Definition: MatrixRepresentationsStatic.h:267
ROOT::Math::MatRepStd::kSize
@ kSize
return no of elements: rows*columns
Definition: MatrixRepresentationsStatic.h:109
f
#define f(i)
Definition: RSha256.hxx:104
ROOT::Math::MatRepSym::off
static int off(int i)
Definition: MatrixRepresentationsStatic.h:307
ROOT::Math::MatRepSym
MatRepSym Matrix storage representation for a symmetric matrix of dimension NxN This class is a templ...
Definition: MatrixRepresentationsStatic.h:213
ROOT::Math::rowOffsetsUtils::make_indices
Definition: MatrixRepresentationsStatic.h:171
ROOT::Math::MatRepStd::Array
T * Array()
Definition: MatrixRepresentationsStatic.h:72
F
#define F(x, y, z)
ROOT::Math::MatRepStd::apply
T apply(unsigned int i) const
Definition: MatrixRepresentationsStatic.h:70
ROOT::Math::MatRepStd::operator()
const T & operator()(unsigned int i, unsigned int j) const
Definition: MatrixRepresentationsStatic.h:60
ROOT::Math::MatRepStd::kCols
@ kCols
return no. of matrix columns
Definition: MatrixRepresentationsStatic.h:107
ROOT::Math::MatRepStd::Array
const T * Array() const
Definition: MatrixRepresentationsStatic.h:74
ROOT::Math::MatRepSym::Array
T * Array()
Definition: MatrixRepresentationsStatic.h:240
ROOT::Math::MatRepSym::operator-=
MatRepSym< T, D > & operator-=(const R &)
self subtraction : only symmetric to symmetric allowed
Definition: MatrixRepresentationsStatic.h:276
N
#define N
ROOT::Math::MatRepSym::operator+=
MatRepSym< T, D > & operator+=(const R &)
self addition : only symmetric to symmetric allowed
Definition: MatrixRepresentationsStatic.h:262
ROOT::Math::MatRepStd::operator-=
MatRepStd< T, D1, D2 > & operator-=(const R &rhs)
Definition: MatrixRepresentationsStatic.h:83
v
@ v
Definition: rootcling_impl.cxx:3664
ROOT::Math::MatRepSym::off2
static constexpr int off2(int i, int j)
Definition: MatrixRepresentationsStatic.h:304
ROOT::Math::rowOffsetsUtils::do_make
constexpr std::array< decltype(std::declval< F >)(std::declval< int >))), sizeof...(I)> do_make(F f, indices< I... >)
Definition: MatrixRepresentationsStatic.h:178
ROOT::Math::MatRepSym::operator-=
MatRepSym< T, D > & operator-=(const MatRepSym &rhs)
Definition: MatrixRepresentationsStatic.h:281
ROOT::Math::RowOffsets::fOff
int fOff[D *D]
Definition: MatrixRepresentationsStatic.h:146
ROOT::Math::MatRepStd
Expression wrapper class for Matrix objects.
Definition: MatrixRepresentationsStatic.h:54
ROOT::Math::MatRepSym::operator==
bool operator==(const R &rhs) const
Definition: MatrixRepresentationsStatic.h:286
ROOT::Math::MatRepStd::operator=
MatRepStd< T, D1, D2 > & operator=(const R &rhs)
Definition: MatrixRepresentationsStatic.h:89
R
#define R(a, b, c, d, e, f, g, h, i)
Definition: RSha256.hxx:110
ROOT::Math::RowOffsets
Static structure to keep the conversion from (i,j) to offsets in the storage data for a symmetric mat...
Definition: MatrixRepresentationsStatic.h:131
ROOT::Math::MatRepSym::apply
T apply(unsigned int i) const
Definition: MatrixRepresentationsStatic.h:236
ROOT::Math::MatRepSym::operator=
MatRepSym< T, D > & operator=(const R &)
assignment : only symmetric to symmetric allowed
Definition: MatrixRepresentationsStatic.h:248
ROOT::Math::rowOffsetsUtils::make_indices_impl< N, indices< Indices... >, N >::type
indices< Indices... > type
Definition: MatrixRepresentationsStatic.h:167
ROOT::Math::RowOffsets::RowOffsets
RowOffsets()
Definition: MatrixRepresentationsStatic.h:132
ROOT::Math::MatRepSym::kCols
@ kCols
return no. of matrix columns
Definition: MatrixRepresentationsStatic.h:298
ROOT::Math::MatRepSym::Array
const T * Array() const
Definition: MatrixRepresentationsStatic.h:242
ROOT::Math::MatRepSym::fArray
T fArray[kSize]
Definition: MatrixRepresentationsStatic.h:322
ROOT::Math::MatRepSym::operator[]
T const & operator[](unsigned int i) const
Definition: MatrixRepresentationsStatic.h:232
ROOT::Math::RowOffsets::apply
int apply(unsigned int i) const
Definition: MatrixRepresentationsStatic.h:145
ROOT::Math::MatRepSym::operator()
T const & operator()(unsigned int i, unsigned int j) const
Definition: MatrixRepresentationsStatic.h:225
ROOT::Math::MatRepSym::value_type
T value_type
Definition: MatrixRepresentationsStatic.h:219
ROOT::Math::RowOffsets::operator()
int operator()(unsigned int i, unsigned int j) const
Definition: MatrixRepresentationsStatic.h:144
ROOT::Math::MatRepSym::offset
static constexpr unsigned int offset(unsigned int i, unsigned int j)
Definition: MatrixRepresentationsStatic.h:313
ROOT::Math::MatRepStd::operator==
bool operator==(const R &rhs) const
Definition: MatrixRepresentationsStatic.h:95
ROOT::Math::MatRepSym::operator[]
T & operator[](unsigned int i)
Definition: MatrixRepresentationsStatic.h:228
ROOT::Math::MatRepSym::kSize
@ kSize
return no of elements: rows*columns
Definition: MatrixRepresentationsStatic.h:300
ROOT::Math::rowOffsetsUtils::indices
Definition: MatrixRepresentationsStatic.h:153
ROOT::Math::MatRepSym::off1
static constexpr int off1(int i)
Definition: MatrixRepresentationsStatic.h:305
ROOT::Math::MatRepStd::kRows
@ kRows
return no. of matrix rows
Definition: MatrixRepresentationsStatic.h:105
ROOT::Math::MatRepSym::operator=
MatRepSym< T, D > & operator=(const MatRepSym &rhs)
Definition: MatrixRepresentationsStatic.h:253
StaticCheck.h
ROOT::Math::MatRepStd::value_type
T value_type
Definition: MatrixRepresentationsStatic.h:58
ROOT::Math::rowOffsetsUtils::make_indices_impl< I, indices< Indices... >, N >::type
make_indices_impl< I+1, indices< Indices..., I >, N >::type type
Definition: MatrixRepresentationsStatic.h:162
ROOT::Math::MatRepStd::fArray
T fArray[kSize]
Definition: MatrixRepresentationsStatic.h:114
ROOT::Math::MatRepStd::operator[]
T & operator[](unsigned int i)
Definition: MatrixRepresentationsStatic.h:66
ROOT::Math::Chebyshev::T
double T(double x)
Definition: ChebyshevPol.h:34
ROOT::Math::MatRepSym::operator()
T & operator()(unsigned int i, unsigned int j)
Definition: MatrixRepresentationsStatic.h:222
I
#define I(x, y, z)
STATIC_CHECK
#define STATIC_CHECK(expr, msg)
Definition: StaticCheck.h:56
ROOT::Math::MatRepSym::kRows
@ kRows
return no. of matrix rows
Definition: MatrixRepresentationsStatic.h:296
ROOT
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...
Definition: EExecutionPolicy.hxx:4
Math
Namespace for new Math classes and functions.