Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
StandardMatrix.hpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File: StandardMatrix.hpp
4 //
5 // For more information, please see: http://www.nektar.info
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: Interface classes for matrices
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #ifndef NEKTAR_LIB_UTILITIES_LINEAR_ALGEBRA_STANDARD_MATRIX_HPP
37 #define NEKTAR_LIB_UTILITIES_LINEAR_ALGEBRA_STANDARD_MATRIX_HPP
38 
39 #include <LibUtilities/LibUtilitiesDeclspec.h>
40 
41 #include <LibUtilities/LinearAlgebra/MatrixBase.hpp>
42 #include <LibUtilities/LinearAlgebra/PointerWrapper.h>
43 #include <LibUtilities/LinearAlgebra/NekVectorFwd.hpp>
44 #include <ExpressionTemplates/ExpressionTemplates.hpp>
45 #include <LibUtilities/LinearAlgebra/NekMatrixMetadata.hpp>
46 #include <LibUtilities/LinearAlgebra/MatrixFuncs.h>
47 #include <LibUtilities/LinearAlgebra/MatrixSize.hpp>
48 
49 #include <LibUtilities/BasicUtils/Vmath.hpp>
50 #include <LibUtilities/BasicUtils/VmathArray.hpp>
51 
52 #include <ExpressionTemplates/Node.hpp>
53 #include <ExpressionTemplates/ExpressionEvaluator.hpp>
54 
55 namespace Nektar
56 {
57  /// \brief Standard Matrix
58  /// \param DataType The type stored in each element.
59  ///
60  /// Matrices are stored in column major order to make it easier to interoperate with
61  /// Blas and Lapack.
62  template<typename DataType>
63  class NekMatrix<DataType, StandardMatrixTag> : public Matrix<DataType>
64  {
65  public:
66  typedef Matrix<DataType> BaseType;
67  typedef NekMatrix<DataType, StandardMatrixTag> ThisType;
68  typedef DataType NumberType;
69  typedef const DataType& ConstGetValueType;
70  typedef DataType GetValueType;
71 
72  public:
73 
74  template<typename T, typename MatrixType>
75  class iterator_impl
76  {
77  public:
78  typedef T value_type;
79 
80  template<typename Z>
81  struct TagType
82  {
83  typedef std::forward_iterator_tag type;
84  };
85 
86  template<typename Z>
87  struct TagType<const Z>
88  {
89  typedef std::input_iterator_tag type;
90  };
91 
92  typedef typename TagType<T>::type iterator_category;
93  typedef unsigned int difference_type;
94  typedef typename boost::call_traits<value_type>::reference reference;
95  typedef typename boost::call_traits<value_type>::const_reference const_reference;
96  typedef typename boost::remove_reference<value_type>::type* pointer;
97 
98  public:
99  iterator_impl(pointer d, pointer e, bool isEnd = false) :
100  m_data(d),
101  m_end(e),
102  m_curRow(std::numeric_limits<unsigned int>::max()),
103  m_curColumn(std::numeric_limits<unsigned int>::max()),
104  m_matrix(NULL),
105  m_curIndex(std::numeric_limits<unsigned int>::max()),
106  m_transpose('N')
107  {
108  if( isEnd )
109  {
110  m_data = m_end;
111  }
112  }
113 
114  iterator_impl(MatrixType* m, char transpose, bool isEnd = false) :
115  m_data(NULL),
116  m_end(NULL),
117  m_curRow(0),
118  m_curColumn(0),
119  m_matrix(m),
120  m_curIndex(0),
121  m_transpose(transpose)
122  {
123  if( isEnd )
124  {
125  m_curRow = std::numeric_limits<unsigned int>::max();
126  m_curColumn = std::numeric_limits<unsigned int>::max();
127  m_curIndex = std::numeric_limits<unsigned int>::max();
128  }
129  }
130 
131  iterator_impl(const iterator_impl<T, MatrixType>& rhs) :
132  m_data(rhs.m_data),
133  m_end(rhs.m_end),
134  m_curRow(rhs.m_curRow),
135  m_curColumn(rhs.m_curColumn),
136  m_matrix(rhs.m_matrix),
137  m_curIndex(rhs.m_curIndex),
138  m_transpose(rhs.m_transpose)
139  {
140  }
141 
142  iterator_impl<T, MatrixType>& operator=(const iterator_impl<T, MatrixType>& rhs)
143  {
144  m_data = rhs.m_data;
145  m_end = rhs.m_end;
146  m_curRow = rhs.m_curRow;
147  m_curColumn = rhs.m_curColumn;
148  m_matrix = rhs.m_matrix;
149  m_curIndex = rhs.m_curIndex;
150  m_transpose = rhs.m_transpose;
151  return *this;
152  }
153 
154  reference operator*()
155  {
156  if( m_data )
157  {
158  ASSERTL1(m_data < m_end, "Attempt to dereference matrix iterator after its end.");
159  return *m_data;
160  }
161  else
162  {
163  return m_matrix->GetPtr()[m_curIndex];
164  }
165  }
166 
167  const_reference operator*() const
168  {
169  if( m_data )
170  {
171  ASSERTL1(m_data < m_end, "Attempt to dereference matrix iterator after its end.");
172  return *m_data;
173  }
174  else
175  {
176  return m_matrix->GetPtr()[m_curIndex];
177  }
178  }
179 
180  /// \brief Prefix increment operator.
181  iterator_impl<T, MatrixType>& operator++()
182  {
183  if( m_data )
184  {
185  ++m_data;
186  }
187  else
188  {
189  boost::tie(m_curRow, m_curColumn) =
190  m_matrix->Advance(m_curRow, m_curColumn, m_transpose);
191  if( m_curRow == std::numeric_limits<unsigned int>::max() )
192  {
193  m_curIndex = m_curRow;
194  }
195  else
196  {
197  m_curIndex = m_matrix->CalculateIndex(m_curRow, m_curColumn, m_transpose);
198  }
199  }
200  return *this;
201  }
202 
203  /// \postfix increment operator.
204  iterator_impl<T, MatrixType> operator++(int)
205  {
206  iterator_impl<T, MatrixType> result = *this;
207  ++(*this);
208  return result;
209  }
210 
211  bool operator==(const iterator_impl<T, MatrixType>& rhs)
212  {
213  return m_data == rhs.m_data &&
214  m_end == rhs.m_end &&
215  m_curRow == rhs.m_curRow &&
216  m_curColumn == rhs.m_curColumn &&
217  m_matrix == rhs.m_matrix &&
218  m_curIndex == rhs.m_curIndex &&
219  m_transpose == rhs.m_transpose;
220  }
221 
222  bool operator!=(const iterator_impl<T, MatrixType>& rhs)
223  {
224  return !(*this == rhs);
225  }
226 
227  private:
228  // Used when the matrix is not transposed
229  T* m_data;
230  T* m_end;
231 
232  // Used when the matrix is transposed.
233  unsigned int m_curRow;
234  unsigned int m_curColumn;
235  MatrixType* m_matrix;
236  unsigned int m_curIndex;
237  char m_transpose;
238  };
239 
240  public:
241  /// \brief Creates an empty matrix.
242  LIB_UTILITIES_EXPORT NekMatrix() ;
243 
244  /// \brief Creates a rows by columns matrix.
245  /// \brief rows The number of rows in the matrix.
246  /// \brief columns The number of columns in the matrix.
247  LIB_UTILITIES_EXPORT NekMatrix(unsigned int rows, unsigned int columns, MatrixStorage policy = eFULL,
248  unsigned int subDiagonals = std::numeric_limits<unsigned int>::max(),
249  unsigned int superDiagonals = std::numeric_limits<unsigned int>::max());
250 
251  LIB_UTILITIES_EXPORT NekMatrix(unsigned int rows, unsigned int columns,
252  MatrixStorage policy,
253  unsigned int subDiagonals,
254  unsigned int superDiagonals,
255  unsigned int capacity);
256 
257  /// \brief Creates a rows by columns matrix.
258  /// \brief rows The number of rows in the matrix.
259  /// \brief columns The number of columns in the matrix.
260  /// \brief initValue The value used to initialize each element.
261  /// \brief policySpecificData Data the is specific to the storage type assigned to this matrix.
262  /// Check MatrixStoragePolicy<StorageType> to see if the StoragePolicy
263  /// for this matrix has policy specific data.
264  LIB_UTILITIES_EXPORT NekMatrix(unsigned int rows, unsigned int columns, typename boost::call_traits<DataType>::const_reference initValue,
265  MatrixStorage policy = eFULL,
266  unsigned int subDiagonals = std::numeric_limits<unsigned int>::max(),
267  unsigned int superDiagonals = std::numeric_limits<unsigned int>::max());
268 
269  /// \brief Creates a rows by columns matrix.
270  /// \brief rows The number of rows in the matrix.
271  /// \brief columns The number of columns in the matrix.
272  /// \brief data An array of data use to initialize the matrix.
273  /// \brief policySpecificData Data the is specific to the storage type assigned to this matrix.
274  /// Check MatrixStoragePolicy<StorageType> to see if the StoragePolicy
275  /// for this matrix has policy specific data.
276  LIB_UTILITIES_EXPORT NekMatrix(unsigned int rows, unsigned int columns, const DataType* data,
277  MatrixStorage policy = eFULL,
278  unsigned int subDiagonals = std::numeric_limits<unsigned int>::max(),
279  unsigned int superDiagonals = std::numeric_limits<unsigned int>::max());
280 
281  /// \brief Creates a rows by columns matrix.
282  /// \brief rows The number of rows in the matrix.
283  /// \brief columns The number of columns in the matrix.
284  /// \brief d An array of data used to initialize the matrix. Values from d are copied into the matrix.
285  /// \brief policySpecificData Data the is specific to the storage type assigned to this matrix.
286  /// Check MatrixStoragePolicy<StorageType> to see if the StoragePolicy
287  /// for this matrix has policy specific data.
288  LIB_UTILITIES_EXPORT NekMatrix(unsigned int rows, unsigned int columns, const Array<OneD, const DataType>& d,
289  MatrixStorage policy = eFULL,
290  unsigned int subDiagonals = std::numeric_limits<unsigned int>::max(),
291  unsigned int superDiagonals = std::numeric_limits<unsigned int>::max());
292 
293  /// \brief Creates a rows by columns matrix.
294  /// \brief rows The number of rows in the matrix.
295  /// \brief columns The number of columns in the matrix.
296  /// \brief d An array of data used to initialize the matrix. If wrapperType is eCopy, then
297  /// each element is copied from d to the matrix. If wrapperType is eWrapper, then
298  /// the matrix uses d directly as its matrix data and no copies are made.
299  /// \brief policySpecificData Data the is specific to the storage type assigned to this matrix.
300  /// Check MatrixStoragePolicy<StorageType> to see if the StoragePolicy
301  /// for this matrix has policy specific data.
302  LIB_UTILITIES_EXPORT NekMatrix(unsigned int rows, unsigned int columns, const Array<OneD, DataType>& d, PointerWrapper wrapperType = eCopy,
303  MatrixStorage policy = eFULL,
304  unsigned int subDiagonals = std::numeric_limits<unsigned int>::max(),
305  unsigned int superDiagonals = std::numeric_limits<unsigned int>::max());
306 
307  LIB_UTILITIES_EXPORT NekMatrix(const ThisType& rhs);
308 
309  #ifdef NEKTAR_USE_EXPRESSION_TEMPLATES
310  template<typename L, typename Op, typename R>
311  NekMatrix(const expt::Node<L, Op, R>& rhs) :
312  BaseType(0, 0),
313  m_wrapperType(eCopy),
314  m_storagePolicy(eFULL),
315  m_numberOfSuperDiagonals(std::numeric_limits<unsigned int>::max()),
316  m_numberOfSubDiagonals(std::numeric_limits<unsigned int>::max()),
317  m_tempSpace()
318  {
319  boost::tuple<unsigned int, unsigned int, unsigned int> sizes =
320  MatrixSize<expt::Node<L, Op, R>, typename expt::Node<L, Op, R>::Indices, 0>::GetRequiredSize(rhs.GetData());
321 
322  unsigned int rows = sizes.get<0>();
323  unsigned int columns = sizes.get<1>();
324  unsigned int bufferCapacity = sizes.get<2>();
325 
326  m_data = Array<OneD, DataType>(bufferCapacity);
327  this->Resize(rows, columns);
328 
329  expt::ExpressionEvaluator::EvaluateWithoutAliasingCheck(rhs, *this);
330  this->RemoveExcessCapacity();
331  }
332 
333  #endif //NEKTAR_USE_EXPRESSION_TEMPLATES
334 
335 
336  LIB_UTILITIES_EXPORT MatrixStorage GetType() const;
337 
338  LIB_UTILITIES_EXPORT ThisType& operator=(const ThisType& rhs);
339 
340  template<typename InnerMatrixType>
341  ThisType& operator=(const NekMatrix<InnerMatrixType, ScaledMatrixTag>& rhs)
342  {
343  BaseType::operator=(rhs);
344  m_storagePolicy = rhs.GetType();
345  m_numberOfSubDiagonals = rhs.GetNumberOfSubDiagonals();
346  m_numberOfSuperDiagonals = rhs.GetNumberOfSuperDiagonals();
347 
348  ResizeDataArrayIfNeeded();
349 
350  unsigned int requiredStorageSize = GetRequiredStorageSize();
351  DataType scale = rhs.Scale();
352 
353  DataType* lhs_array = m_data.data();
354  const DataType* rhs_array = rhs.GetRawPtr();
355 
356  for(unsigned int i = 0; i < requiredStorageSize; ++i)
357  {
358  lhs_array[i] = scale*rhs_array[i];
359  }
360 
361  return *this;
362 
363  }
364 
365 
366  #ifdef NEKTAR_USE_EXPRESSION_TEMPLATES
367  template<typename L, typename Op, typename R>
368  ThisType& operator=(const expt::Node<L, Op, R>& rhs)
369  {
370  if( this->GetWrapperType() == eWrapper ||
371  expt::ExpressionEvaluator::ContainsAlias(rhs, *this) )
372  {
373  // Assignment into a wrapped matrix can't change the matrix size (the matrix
374  // can be wrapping a portion of a larger array that can't be resized).
375  // If the matrix is wrapped, we'll evaluate the expression into temporary
376  // and then assign the temporary, relying on the operator= for the error
377  // checking relating to size.
378 
379  // If the expression is aliased, we need to do the same thing, as
380  // we can't change the size of the accumulator before evaluation.
381  ThisType temp(rhs);
382  *this = temp;
383  }
384  else
385  {
386  // If the matrix is not wrapped, then we are free to resize as necessary.
387 
388  boost::tuple<unsigned int, unsigned int, unsigned int> sizes =
389  MatrixSize<expt::Node<L, Op, R>, typename expt::Node<L, Op, R>::Indices, 0>::GetRequiredSize(rhs.GetData());
390  unsigned int rows = sizes.get<0>();
391  unsigned int columns = sizes.get<1>();
392  unsigned int bufferSize = sizes.get<2>();
393 
394  if( this->GetRows() != rows ||
395  this->GetColumns() != columns )
396  {
397  this->Resize(rows, columns);
398  }
399 
400  // We don't resize the temp workspace since we don't know if we will
401  // need it. It will be sized when needed during expression evaluation.
402  this->ResizeDataArrayIfNeeded(bufferSize);
403 
404  this->SetTransposeFlag('N');
405 
406  expt::ExpressionEvaluator::Evaluate(rhs, *this);
407  this->RemoveExcessCapacity();
408  }
409 
410  return *this;
411  }
412  #endif //NEKTAR_USE_EXPRESSION_TEMPLATES
413 
414  /// \brief Returns the element value at the given row and column.
415  LIB_UTILITIES_EXPORT ConstGetValueType operator()(unsigned int row, unsigned int column) const;
416 
417  class Proxy
418  {
419  public:
420  Proxy() : m_value(defaultReturnValue) {}
421  explicit Proxy(DataType& value) : m_value(value) {}
422  Proxy(const Proxy& rhs) : m_value(rhs.m_value) {}
423  Proxy& operator=(const Proxy& rhs)
424  {
425  m_value = rhs.m_value;
426  return *this;
427  }
428 
429  DataType& operator*() { return m_value; }
430  const DataType& operator*() const { return m_value; }
431  operator DataType&() { return m_value; }
432  operator const DataType&() const { return m_value; }
433  void operator=(const DataType& newValue)
434  {
435  if( &m_value != &defaultReturnValue )
436  {
437  m_value = newValue;
438  }
439  }
440 
441  private:
442  DataType& m_value;
443  static DataType defaultReturnValue;
444  };
445 
446  /// \brief Returns the element value at the given row and column.
447  /// \brief row The element's row.
448  /// \brief column The element's column.
449  /// \brief transpose If transpose = 'N', then the return value is element [row, column].
450  /// If transpose = 'T', then the return value is element [column, row].
451  LIB_UTILITIES_EXPORT ConstGetValueType operator()(unsigned int row, unsigned int column, char transpose) const;
452 
453  LIB_UTILITIES_EXPORT unsigned int GetRequiredStorageSize() const;
454 
455  LIB_UTILITIES_EXPORT unsigned int CalculateIndex(unsigned int row, unsigned int col, const char transpose) const;
456 
457  /// \brief Returns the element value at the given row and column.
458  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::const_reference GetValue(unsigned int row, unsigned int column) const;
459 
460  /// \brief Returns the element value at the given row and column.
461  /// \brief row The element's row.
462  /// \brief column The element's column.
463  /// \brief transpose If transpose = 'N', then the return value is element [row, column].
464  /// If transpose = 'T', then the return value is element [column, row].
465  LIB_UTILITIES_EXPORT ConstGetValueType GetValue(unsigned int row, unsigned int column, char transpose) const;
466 
467  LIB_UTILITIES_EXPORT const Array<OneD, const DataType>& GetPtr() const;
468 
469  /// \brief Returns the scaling used by this matrix.
470  ///
471  /// Since this matrix is not a scaled matrix, this method always returns 1.
472  LIB_UTILITIES_EXPORT DataType Scale() const;
473 
474  LIB_UTILITIES_EXPORT const DataType* GetRawPtr() const;
475 
476  typedef iterator_impl<const DataType, const ThisType> const_iterator;
477  typedef iterator_impl<DataType, ThisType> iterator;
478 
479  LIB_UTILITIES_EXPORT iterator begin();
480 
481  LIB_UTILITIES_EXPORT iterator begin(char transpose);
482 
483  LIB_UTILITIES_EXPORT iterator end();
484  LIB_UTILITIES_EXPORT iterator end(char transpose);
485 
486  LIB_UTILITIES_EXPORT const_iterator begin() const;
487 
488  LIB_UTILITIES_EXPORT const_iterator begin(char transpose) const;
489 
490  LIB_UTILITIES_EXPORT const_iterator end() const;
491 
492  LIB_UTILITIES_EXPORT const_iterator end(char transpose) const;
493 
494 
495  LIB_UTILITIES_EXPORT unsigned int GetStorageSize() const;
496 
497  // Banded matrices only.
498  // Get the specified number of sub diagonals, or calculate the
499  // number if this object has not been initialized.
500  LIB_UTILITIES_EXPORT unsigned int GetNumberOfSubDiagonals() const;
501 
502  LIB_UTILITIES_EXPORT unsigned int GetNumberOfSuperDiagonals() const;
503 
504  LIB_UTILITIES_EXPORT unsigned int CalculateNumberOfRows() const;
505 
506  /// \brief Returns true if the this matrix and rhs are equivalent.
507  ///
508  /// Two matrices are equivalent if they have the same size and each element
509  /// is the same.
510  LIB_UTILITIES_EXPORT bool operator==(const NekMatrix<DataType, StandardMatrixTag>& rhs) const;
511 
512  LIB_UTILITIES_EXPORT PointerWrapper GetWrapperType() const;
513 
514 
515  LIB_UTILITIES_EXPORT char GetTransposeFlag() const ;
516 
517 
518  LIB_UTILITIES_EXPORT boost::tuples::tuple<unsigned int, unsigned int>
519  Advance(unsigned int curRow, unsigned int curColumn) const;
520 
521  LIB_UTILITIES_EXPORT boost::tuples::tuple<unsigned int, unsigned int>
522  Advance(unsigned int curRow, unsigned int curColumn, char transpose) const;
523 
524  LIB_UTILITIES_EXPORT static ThisType CreateWrapper(const ThisType& rhs);
525 
526  LIB_UTILITIES_EXPORT static boost::shared_ptr<ThisType> CreateWrapper(const boost::shared_ptr<ThisType>& rhs);
527 
528  LIB_UTILITIES_EXPORT void SetSize(unsigned int rows, unsigned int cols);
529 
530  LIB_UTILITIES_EXPORT Proxy operator()(unsigned int row, unsigned int column);
531 
532 
533  LIB_UTILITIES_EXPORT Proxy operator()(unsigned int row, unsigned int column, char transpose);
534 
535  LIB_UTILITIES_EXPORT void SetValue(unsigned int row, unsigned int column, typename boost::call_traits<DataType>::const_reference d);
536 
537  LIB_UTILITIES_EXPORT void SetValue(unsigned int row, unsigned int column, typename boost::call_traits<DataType>::const_reference d, char transpose);
538 
539  LIB_UTILITIES_EXPORT Array<OneD, DataType>& GetPtr();
540 
541  LIB_UTILITIES_EXPORT DataType* GetRawPtr();
542 
543 
544 
545  LIB_UTILITIES_EXPORT NekDouble AbsMaxtoMinEigenValueRatio(void);
546 
547  LIB_UTILITIES_EXPORT void EigenSolve(Array<OneD, NekDouble> &EigValReal,
548  Array<OneD, NekDouble> &EigValImag,
549  Array<OneD, NekDouble> &EigVecs =NullNekDouble1DArray);
550 
551  LIB_UTILITIES_EXPORT void Invert();
552 
553  LIB_UTILITIES_EXPORT Array<OneD, DataType>& GetTempSpace();
554 
555  LIB_UTILITIES_EXPORT void SwapTempAndDataBuffers();
556 
557  LIB_UTILITIES_EXPORT ThisType& operator*=(const NumberType& s);
558 
559  #ifdef NEKTAR_USE_EXPRESSION_TEMPLATES
560  expt::Node<expt::Node<NekMatrix<DataType, StandardMatrixTag> >, expt::NegateOp, void > operator-() const
561  {
562  expt::Node<NekMatrix<DataType, StandardMatrixTag> > leafNode(*this);
563  return expt::Node<expt::Node<NekMatrix<DataType, StandardMatrixTag> >, expt::NegateOp, void >(leafNode);
564  }
565  #else
566  LIB_UTILITIES_EXPORT NekMatrix<DataType, StandardMatrixTag> operator-() const;
567  #endif
568 
569  protected:
570  LIB_UTILITIES_EXPORT NekMatrix(const ThisType& rhs, PointerWrapper wrapperType);
571 
572  LIB_UTILITIES_EXPORT Array<OneD, DataType>& GetData();
573 
574  LIB_UTILITIES_EXPORT void RemoveExcessCapacity();
575 
576  LIB_UTILITIES_EXPORT void ResizeDataArrayIfNeeded(unsigned int requiredStorageSize);
577 
578  LIB_UTILITIES_EXPORT void ResizeDataArrayIfNeeded();
579 
580 
581  private:
582  LIB_UTILITIES_EXPORT void PerformCopyConstruction(const ThisType& rhs);
583 
584  LIB_UTILITIES_EXPORT virtual typename boost::call_traits<DataType>::value_type v_GetValue(unsigned int row, unsigned int column) const;
585 
586  LIB_UTILITIES_EXPORT virtual unsigned int v_GetStorageSize() const ;
587 
588  LIB_UTILITIES_EXPORT virtual MatrixStorage v_GetStorageType() const;
589 
590  // We need to rethink class structure a little. This shouldn't be necessary.
591  LIB_UTILITIES_EXPORT virtual void v_SetValue(unsigned int row, unsigned int column, typename boost::call_traits<DataType>::const_reference d);
592 
593  Array<OneD, DataType> m_data;
594  PointerWrapper m_wrapperType;
595  MatrixStorage m_storagePolicy;
596 
597  // Only used by banded matrices.
598  unsigned int m_numberOfSuperDiagonals;
599  unsigned int m_numberOfSubDiagonals;
600 
601  // Used during chained matrix multiplication as workspace.
602  Array<OneD, DataType> m_tempSpace;
603  };
604 
605 
606  template<typename DataType>
607  DataType NekMatrix<DataType, StandardMatrixTag>::Proxy::defaultReturnValue;
608 
609  template<typename DataType>
610  LIB_UTILITIES_EXPORT NekMatrix<DataType, StandardMatrixTag>
611  Transpose(NekMatrix<DataType, StandardMatrixTag>& rhs);
612 
613  template<typename DataType>
614  LIB_UTILITIES_EXPORT void NegateInPlace(NekMatrix<DataType, StandardMatrixTag>& v);
615 }
616 
617 #ifdef NEKTAR_USE_EXPRESSION_TEMPLATES
618 namespace expt
619 {
620  template<typename DataType>
621  struct IsAlias<Nektar::NekMatrix<DataType, Nektar::StandardMatrixTag>, Nektar::NekMatrix<DataType, Nektar::StandardMatrixTag> >
622  {
623  static bool Apply(const Nektar::NekMatrix<DataType, Nektar::StandardMatrixTag>& lhs, const Nektar::NekMatrix<DataType, Nektar::StandardMatrixTag>& rhs)
624  {
625  return lhs.GetPtr().Overlaps(rhs.GetPtr());
626  }
627  };
628 
629  template<typename DataType, typename NodeType, typename Indices, unsigned int StartIndex>
630  struct CreateFromTree<Nektar::NekMatrix<DataType, Nektar::StandardMatrixTag>, NodeType, Indices, StartIndex>
631  {
632  template<typename ArgVectorType>
633  static Nektar::NekMatrix<DataType, Nektar::StandardMatrixTag> Apply(const ArgVectorType& tree)
634  {
635  boost::tuple<unsigned int, unsigned int, unsigned int> sizes =
636  Nektar::MatrixSize<NodeType, Indices, StartIndex>::GetRequiredSize(tree);
637 
638  unsigned int rows = sizes.get<0>();
639  unsigned int columns = sizes.get<1>();
640  unsigned int bufferSize = sizes.get<2>();
641 
642  return Nektar::NekMatrix<DataType, Nektar::StandardMatrixTag>(rows, columns,
643  Nektar::eFULL, std::numeric_limits<unsigned int>::max(),
644  std::numeric_limits<unsigned int>::max(), bufferSize);
645  }
646  };
647 }
648 #endif
649 
650 
651 #endif //NEKTAR_LIB_UTILITIES_LINEAR_ALGEBRA_STANDARD_MATRIX_HPP
Generated on Sun Mar 15 2015 00:11:34 for Nektar++ by   doxygen 1.8.1.2