Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
NekVector.hpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File: NekVector.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: Generic N-Dimensional Vector.
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #ifndef NEKTAR_LIB_UTILITIES_NEK_VECTOR_HPP
37 #define NEKTAR_LIB_UTILITIES_NEK_VECTOR_HPP
38 
39 #include <ExpressionTemplates/ExpressionTemplates.hpp>
40 #include <LibUtilities/BasicUtils/ErrorUtil.hpp>
41 #include <LibUtilities/BasicUtils/SharedArray.hpp>
42 
43 #include <LibUtilities/LinearAlgebra/NekPoint.hpp>
44 
45 #include <LibUtilities/LinearAlgebra/NekVectorMetadata.hpp>
46 #include <LibUtilities/LinearAlgebra/NekVectorFwd.hpp>
47 #include <LibUtilities/LinearAlgebra/NekVectorConstantSized.hpp>
48 #include <LibUtilities/LinearAlgebra/NekVectorVariableSized.hpp>
49 #include <LibUtilities/LinearAlgebra/PointerWrapper.h>
50 
51 #include <LibUtilities/LinearAlgebra/MatrixSize.hpp>
52 
53 #include <LibUtilities/BasicUtils/OperatorGenerators.hpp>
54 
55 #include <functional>
56 #include <algorithm>
57 #include <math.h>
58 
59 #include <boost/call_traits.hpp>
60 #include <boost/type_traits.hpp>
61 #include <boost/shared_array.hpp>
62 #include <boost/typeof/typeof.hpp>
63 
64 #include BOOST_TYPEOF_INCREMENT_REGISTRATION_GROUP()
65 
66 namespace Nektar
67 {
68  template<typename DataType>
69  class NekVector
70  {
71  public:
72  /// \brief Creates an empty vector.
73  LIB_UTILITIES_EXPORT NekVector();
74 
75  /// \brief Creates a vector of given size. The elements are not initialized.
76  LIB_UTILITIES_EXPORT explicit NekVector(unsigned int size);
77 
78  /// \brief Creates a vector with given size and initial value.
79  LIB_UTILITIES_EXPORT NekVector(unsigned int size, typename boost::call_traits<DataType>::const_reference a);
80 
81 
82  LIB_UTILITIES_EXPORT explicit NekVector(const std::string& vectorValues);
83 
84  LIB_UTILITIES_EXPORT NekVector(typename boost::call_traits<DataType>::const_reference x,
85  typename boost::call_traits<DataType>::const_reference y,
86  typename boost::call_traits<DataType>::const_reference z);
87 
88  LIB_UTILITIES_EXPORT NekVector(const NekVector<DataType>& rhs);
89 
90  LIB_UTILITIES_EXPORT NekVector(unsigned int size, const DataType* const ptr);
91  LIB_UTILITIES_EXPORT explicit NekVector(const Array<OneD, DataType>& ptr, PointerWrapper h = eCopy);
92  LIB_UTILITIES_EXPORT NekVector(unsigned int size, Array<OneD, DataType>& ptr, PointerWrapper h = eCopy);
93 
94  LIB_UTILITIES_EXPORT NekVector(unsigned int size, const Array<OneD, const DataType>& ptr, PointerWrapper h = eCopy);
95 
96  #ifdef NEKTAR_USE_EXPRESSION_TEMPLATES
97  template<typename L, typename Op, typename R>
98  NekVector(const expt::Node<L, Op, R>& rhs) :
99  m_size(MatrixSize<expt::Node<L, Op, R>, typename expt::Node<L, Op, R>::Indices, 0>::template GetRequiredSize(rhs.GetData()).template get<0>()),
100  m_data(m_size),
101  m_wrapperType(eCopy)
102  {
103  expt::ExpressionEvaluator::EvaluateWithoutAliasingCheck(rhs, *this);
104  }
105 
106  template<typename L, typename Op, typename R>
107  NekVector<DataType>& operator=(const expt::Node<L, Op, R>& rhs)
108  {
109  boost::tuple<unsigned int, unsigned int, unsigned int> sizes =
110  MatrixSize<expt::Node<L, Op, R>, typename expt::Node<L, Op, R>::Indices, 0>::GetRequiredSize(rhs.GetData());
111  unsigned int newRows = sizes.get<0>();
112 
113  this->SetSize(newRows);
114  if( this->GetWrapperType() == eCopy )
115  {
116  if( this->GetData().num_elements() < newRows )
117  {
118  this->SetData(Array<OneD, DataType>(newRows));
119  }
120  }
121  else
122  {
123  ASSERTL0(this->GetData().num_elements() >= newRows,
124  "Attempting to store too many elements in a wrapped vector.");
125  }
126 
127  expt::ExpressionEvaluator::Evaluate(rhs, *this);
128  return *this;
129  }
130  #endif //NEKTAR_USE_EXPRESSION_TEMPLATES
131 
132  LIB_UTILITIES_EXPORT ~NekVector();
133 
134  LIB_UTILITIES_EXPORT NekVector<DataType>& operator=(const NekVector<DataType>& rhs);
135 
136 
137  /// \brief Returns the number of dimensions for the point.
138  LIB_UTILITIES_EXPORT unsigned int GetDimension() const;
139 
140  LIB_UTILITIES_EXPORT unsigned int GetRows() const;
141 
142  LIB_UTILITIES_EXPORT DataType* GetRawPtr();
143 
144  LIB_UTILITIES_EXPORT Array<OneD, DataType>& GetPtr();
145 
146  LIB_UTILITIES_EXPORT const DataType* GetRawPtr() const;
147 
148  LIB_UTILITIES_EXPORT const Array<OneD, const DataType>& GetPtr() const;
149 
150  typedef DataType* iterator;
151  LIB_UTILITIES_EXPORT iterator begin();
152  LIB_UTILITIES_EXPORT iterator end();
153 
154  typedef const DataType* const_iterator;
155  LIB_UTILITIES_EXPORT const_iterator begin() const;
156  LIB_UTILITIES_EXPORT const_iterator end() const;
157 
158  /// \brief Returns i^{th} element.
159  /// \param i The element to return.
160  /// \pre i < dim
161  /// \return A reference to the i^{th} element.
162  ///
163  /// Retrieves the i^{th} element. Since it returns a reference you may
164  /// assign a new value (i.e., p(2) = 3.2;)
165  ///
166  /// This operator performs range checking.
167  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::reference operator()(unsigned int i);
168 
169  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::reference operator[](unsigned int i);
170 
171  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::reference x();
172 
173  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::reference y();
174 
175  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::reference z();
176 
177  LIB_UTILITIES_EXPORT void SetX(typename boost::call_traits<DataType>::const_reference val);
178 
179  LIB_UTILITIES_EXPORT void SetY(typename boost::call_traits<DataType>::const_reference val);
180 
181  LIB_UTILITIES_EXPORT void SetZ(typename boost::call_traits<DataType>::const_reference val);
182 
183  LIB_UTILITIES_EXPORT NekVector<DataType>& operator+=(const NekVector<DataType>& rhs);
184 
185  LIB_UTILITIES_EXPORT NekVector<DataType>& operator-=(const NekVector<DataType>& rhs);
186 
187  LIB_UTILITIES_EXPORT NekVector<DataType>& operator*=(typename boost::call_traits<DataType>::const_reference rhs);
188 
189  LIB_UTILITIES_EXPORT NekVector<DataType>& operator/=(typename boost::call_traits<DataType>::const_reference rhs);
190 
191  LIB_UTILITIES_EXPORT void Normalize();
192 
193  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::const_reference operator()(unsigned int i) const;
194 
195  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::const_reference operator[](unsigned int i) const;
196 
197  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::const_reference x() const;
198 
199  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::const_reference y() const;
200 
201  LIB_UTILITIES_EXPORT typename boost::call_traits<DataType>::const_reference z() const;
202 
203  #ifdef NEKTAR_USE_EXPRESSION_TEMPLATES
204  expt::Node<expt::Node<NekVector<DataType> >, expt::NegateOp, void > operator-() const
205  {
206  expt::Node<NekVector<DataType> > leafNode(*this);
207  return expt::Node<expt::Node<NekVector<DataType> >, expt::NegateOp, void >(leafNode);
208  }
209  #else
210  LIB_UTILITIES_EXPORT NekVector<DataType> operator-() const;
211  #endif
212 
213  LIB_UTILITIES_EXPORT DataType Magnitude() const;
214  LIB_UTILITIES_EXPORT DataType Dot(const NekVector<DataType>& rhs) const;
215 
216  LIB_UTILITIES_EXPORT NekVector<DataType> Cross(const NekVector<DataType>& rhs) const;
217 
218  LIB_UTILITIES_EXPORT std::string AsString() const;
219 
220  // Norms
221  LIB_UTILITIES_EXPORT DataType L1Norm() const;
222  LIB_UTILITIES_EXPORT DataType L2Norm() const;
223  LIB_UTILITIES_EXPORT DataType InfinityNorm() const;
224 
225 
226  LIB_UTILITIES_EXPORT PointerWrapper GetWrapperType() const;
227 
228  protected:
229 
230  LIB_UTILITIES_EXPORT Array<OneD, DataType>& GetData();
231  LIB_UTILITIES_EXPORT void SetSize(unsigned int s);
232  LIB_UTILITIES_EXPORT void SetWrapperType(PointerWrapper p);
233  LIB_UTILITIES_EXPORT void SetData(const Array<OneD, DataType>& newData);
234  LIB_UTILITIES_EXPORT void Resize(unsigned int newSize);
235 
236  private:
237  // Prevents accidental use of wrapped mode around ConstArrays.
238  NekVector(const Array<OneD, const DataType>& ptr, PointerWrapper h);
239 
240  unsigned int m_size;
241  Array<OneD, DataType> m_data;
242  PointerWrapper m_wrapperType;
243  };
244 
245  BOOST_TYPEOF_REGISTER_TEMPLATE(NekVector, 1);
246 
247  template<typename DataType>
248  LIB_UTILITIES_EXPORT void Add(NekVector<DataType>& result,
249  const NekVector<DataType>& lhs,
250  const NekVector<DataType>& rhs);
251 
252  template<typename DataType>
253  LIB_UTILITIES_EXPORT void AddNegatedLhs(NekVector<DataType>& result,
254  const NekVector<DataType>& lhs,
255  const NekVector<DataType>& rhs);
256 
257  template<typename DataType>
258  LIB_UTILITIES_EXPORT void AddEqual(NekVector<DataType>& result,
259  const NekVector<DataType>& rhs);
260 
261  template<typename DataType>
262  LIB_UTILITIES_EXPORT void AddEqualNegatedLhs(NekVector<DataType>& result,
263  const NekVector<DataType>& rhs);
264 
265  template<typename LhsDataType,
266  typename RhsDataType>
267  LIB_UTILITIES_EXPORT NekVector<LhsDataType> Add(const NekVector<LhsDataType>& lhs,
268  const NekVector<RhsDataType>& rhs);
269 
270 
271 
272  template<typename ResultDataType, typename InputDataType>
273  LIB_UTILITIES_EXPORT void Subtract(NekVector<ResultDataType>& result,
274  const NekVector<InputDataType>& lhs,
275  const NekVector<InputDataType>& rhs);
276 
277  template<typename ResultDataType, typename InputDataType>
278  LIB_UTILITIES_EXPORT void SubtractNegatedLhs(NekVector<ResultDataType>& result,
279  const NekVector<InputDataType>& lhs,
280  const NekVector<InputDataType>& rhs);
281 
282  template<typename ResultDataType, typename InputDataType>
283  LIB_UTILITIES_EXPORT void SubtractEqual(NekVector<ResultDataType>& result,
284  const NekVector<InputDataType>& rhs);
285 
286  template<typename ResultDataType, typename InputDataType>
287  LIB_UTILITIES_EXPORT void SubtractEqualNegatedLhs(NekVector<ResultDataType>& result,
288  const NekVector<InputDataType>& rhs);
289 
290  template<typename DataType>
291  NekVector<DataType>
292  LIB_UTILITIES_EXPORT Subtract(const NekVector<DataType>& lhs,
293  const NekVector<DataType>& rhs);
294 
295 
296 
297 
298 
299  template<typename ResultDataType, typename InputDataType>
300  void LIB_UTILITIES_EXPORT Divide(NekVector<ResultDataType>& result,
301  const NekVector<InputDataType>& lhs,
302  const NekDouble& rhs);
303 
304  template<typename ResultDataType>
305  void LIB_UTILITIES_EXPORT DivideEqual(NekVector<ResultDataType>& result,
306  const NekDouble& rhs);
307 
308  template<typename DataType>
309  NekVector<DataType>
310  LIB_UTILITIES_EXPORT Divide(const NekVector<DataType>& lhs,
311  const NekDouble& rhs);
312 
313 
314  template<typename ResultDataType, typename InputDataType>
315  void LIB_UTILITIES_EXPORT Multiply(NekVector<ResultDataType>& result,
316  const NekVector<InputDataType>& lhs,
317  const NekVector<InputDataType>& rhs);
318 
319  template<typename ResultDataType, typename InputDataType>
320  void LIB_UTILITIES_EXPORT MultiplyEqual(NekVector<ResultDataType>& result,
321  const NekVector<InputDataType>& rhs);
322 
323  template<typename DataType, typename InputDataType>
324  NekVector<DataType>
325  LIB_UTILITIES_EXPORT Multiply(const NekVector<DataType>& lhs,
326  const NekVector<InputDataType>& rhs);
327 
328  template<typename ResultDataType, typename InputDataType>
329  void LIB_UTILITIES_EXPORT Multiply(NekVector<ResultDataType>& result,
330  const NekVector<InputDataType>& lhs,
331  const NekDouble& rhs);
332 
333  template<typename ResultDataType>
334  void LIB_UTILITIES_EXPORT MultiplyEqual(NekVector<ResultDataType>& result,
335  const NekDouble& rhs);
336 
337  template<typename DataType>
338  NekVector<DataType>
339  LIB_UTILITIES_EXPORT Multiply(const NekVector<DataType>& lhs,
340  const NekDouble& rhs);
341 
342  template<typename ResultDataType, typename InputDataType>
343  void LIB_UTILITIES_EXPORT Multiply(NekVector<ResultDataType>& result,
344  const NekDouble& lhs,
345  const NekVector<InputDataType>& rhs);
346 
347  template<typename ResultDataType, typename InputDataType>
348  void LIB_UTILITIES_EXPORT MultiplyInvertedLhs(NekVector<ResultDataType>& result,
349  const NekDouble& lhs,
350  const NekVector<InputDataType>& rhs);
351 
352  template<typename DataType>
353  NekVector<DataType>
354  LIB_UTILITIES_EXPORT Multiply(const DataType& lhs,
355  const NekVector<DataType>& rhs);
356 
357  GENERATE_MULTIPLICATION_OPERATOR(NekVector, 1, NekDouble, 0);
358  GENERATE_MULTIPLICATION_OPERATOR(NekDouble, 0, NekVector, 1);
359  GENERATE_MULTIPLICATION_OPERATOR(NekVector, 1, NekVector, 1);
360 
361  GENERATE_DIVISION_OPERATOR(NekVector, 1, NekDouble, 0);
362  GENERATE_ADDITION_OPERATOR(NekVector, 1, NekVector, 1);
363  GENERATE_SUBTRACTION_OPERATOR(NekVector, 1, NekVector, 1);
364 
365 
366  template<typename DataType>
367  LIB_UTILITIES_EXPORT std::ostream& operator<<(std::ostream& os, const NekVector<DataType>& rhs);
368 
369  template<typename DataType>
370  LIB_UTILITIES_EXPORT NekVector<DataType> createVectorFromPoints(const NekPoint<DataType>& source, const NekPoint<DataType>& dest);
371 
372  template<typename DataType>
373  LIB_UTILITIES_EXPORT NekPoint<DataType> findPointAlongVector(const NekVector<DataType>& lhs, const DataType& t);
374 
375  template<typename DataType>
376  LIB_UTILITIES_EXPORT bool operator==(const NekVector<DataType>& lhs, const NekVector<DataType>& rhs);
377 
378  template<typename DataType>
379  LIB_UTILITIES_EXPORT bool operator!=(const NekVector<DataType>& lhs, const NekVector<DataType>& rhs);
380 
381  template<typename DataType>
382  LIB_UTILITIES_EXPORT DataType Magnitude(const NekVector<DataType>& v);
383 
384  template<typename DataType>
385  LIB_UTILITIES_EXPORT DataType Dot(const NekVector<DataType>& lhs,
386  const NekVector<DataType>& rhs);
387 
388  template<typename DataType>
389  LIB_UTILITIES_EXPORT std::vector<DataType> FromString(const std::string& str);
390 
391  /// \todo Do the Norms with Blas where applicable.
392  template<typename DataType>
393  LIB_UTILITIES_EXPORT DataType L1Norm(const NekVector<DataType>& v);
394 
395  template<typename DataType>
396  LIB_UTILITIES_EXPORT DataType L2Norm(const NekVector<DataType>& v);
397 
398  template<typename DataType>
399  LIB_UTILITIES_EXPORT DataType InfinityNorm(const NekVector<DataType>& v);
400 
401  template<typename DataType>
402  LIB_UTILITIES_EXPORT NekVector<DataType> Negate(const NekVector<DataType>& v);
403 
404  template<typename DataType>
405  LIB_UTILITIES_EXPORT void NegateInPlace(NekVector<DataType>& v);
406 
407  LIB_UTILITIES_EXPORT void NegateInPlace(NekDouble& v);
408  LIB_UTILITIES_EXPORT void InvertInPlace(NekDouble& v);
409 
410  template<typename DataType>
411  LIB_UTILITIES_EXPORT void Normalize(NekVector<DataType>& v);
412 
413  template<typename DataType>
414  LIB_UTILITIES_EXPORT NekVector<DataType> Cross(const NekVector<DataType>& lhs,
415  const NekVector<DataType>& rhs);
416  template<typename DataType>
417  LIB_UTILITIES_EXPORT std::string AsString(const NekVector<DataType>& v);
418 
419 }
420 
421 #ifdef NEKTAR_USE_EXPRESSION_TEMPLATES
422 namespace expt
423 {
424  template<typename DataType, typename L, typename Op, typename R>
425  DataType Dot(const Nektar::NekVector<DataType>& lhs,
426  const expt::Node<L, Op, R>& expr)
427  {
428  // Evalaute the expression first, expression templates don't chain past
429  // this point since the return value is a scalar.
430  typename expt::Node<L, Op, R>::ResultType rhs = expt::ExpressionEvaluator::Evaluate(expr);
431  return Dot(lhs, rhs);
432  }
433 
434  template<typename DataType>
435  struct IsAlias<Nektar::NekVector<DataType>, Nektar::NekVector<DataType> >
436  {
437  static bool Apply(const Nektar::NekVector<DataType>& lhs, const Nektar::NekVector<DataType>& rhs)
438  {
439  return lhs.GetPtr().Overlaps(rhs.GetPtr());
440  }
441  };
442 
443  template<typename DataType, typename NodeType, typename Indices, unsigned int StartIndex>
444  struct CreateFromTree<Nektar::NekVector<DataType>, NodeType, Indices, StartIndex>
445  {
446  template<typename ArgVectorType>
447  static Nektar::NekVector<DataType> Apply(const ArgVectorType& tree)
448  {
449  boost::tuple<unsigned int, unsigned int, unsigned int> sizes =
450  Nektar::MatrixSize<NodeType, Indices, StartIndex>::GetRequiredSize(tree);
451 
452  unsigned int rows = sizes.get<0>();
453  return Nektar::NekVector<DataType>(rows);
454  }
455  };
456 
457  // Override default expression handling for addition/subtraction of vectors.
458  // Optimal execution is obtained by loop unrolling.
459 
460  template<typename NodeType, typename enabled = void>
461  struct NodeCanUnroll : public boost::false_type {};
462 
463  template<typename Type>
464  struct NodeCanUnroll<expt::Node<Type, void, void>,
465  typename boost::enable_if
466  <
467  Nektar::IsVector<typename expt::Node<Type, void, void>::ResultType>
468  >::type > : public boost::true_type
469  {
470  };
471 
472  template<typename LhsType, typename OpType, typename RhsType>
473  struct NodeCanUnroll<expt::Node<LhsType, OpType, RhsType>,
474  typename boost::enable_if
475  <
476  boost::mpl::and_
477  <
478  Nektar::IsVector<typename LhsType::ResultType>,
479  Nektar::IsVector<typename RhsType::ResultType>,
480  NodeCanUnroll<LhsType>,
481  NodeCanUnroll<RhsType>,
482  //boost::mpl::or_
483  //<
484  boost::is_same<OpType, expt::AddOp>
485  // boost::is_same<OpType, expt::SubtractOp>
486  //>
487  > >::type >: public boost::true_type
488  {
489  };
490 
491  template<typename NodeType, typename IndicesType, unsigned int index>
492  struct Accumulate;
493 
494  template<typename LhsType, typename IndicesType, unsigned int index>
495  struct Accumulate<expt::Node<LhsType, void, void>, IndicesType, index>
496  {
497  static const unsigned int MappedIndex = boost::mpl::at_c<IndicesType, index>::type::value;
498 
499  template<typename ResultType, typename ArgumentVectorType>
500  static void Execute(ResultType& accumulator, const ArgumentVectorType& args, unsigned int i)
501  {
502  accumulator = boost::fusion::at_c<MappedIndex>(args)[i];
503  }
504  };
505 
506  template<typename LhsType, typename Op, typename RhsType, typename IndicesType, unsigned int index>
507  struct Accumulate<expt::Node<LhsType, Op, RhsType>, IndicesType, index>
508  {
509  static const int rhsNodeIndex = index + LhsType::TotalCount;
510 
511  template<typename ResultType, typename ArgumentVectorType>
512  static void Execute(ResultType& accumulator, const ArgumentVectorType& args, unsigned int i)
513  {
514  Accumulate<LhsType, IndicesType, index>::Execute(accumulator, args, i);
515  ResultType rhs;
516  Accumulate<RhsType, IndicesType, rhsNodeIndex>::Execute(rhs, args, i);
517  Op::OpEqual(accumulator, rhs);
518  }
519  };
520 
521 
522 
523  template<typename IndicesType, unsigned int startIndex, unsigned int endIndex, typename enabled=void>
524  struct Unroll;
525 
526  #ifndef NEKTAR_NEKVECTOR_MAX_UNROLL_ARGS
527  #define NEKTAR_NEKVECTOR_MAX_UNROLL_ARGS 10
528  #endif
529 
530  #define NEKTAR_NEKVECTOR_UNROLL_GENERATE_INDEX(z, n, IndexName) \
531  static const unsigned int BOOST_PP_CAT(IndexName, n) = boost::mpl::at_c<IndicesType, startIndex+n>::type::value;
532 
533  #define NEKTAR_NEKVECTOR_UNROLL_GENERATE_VARIABLE(z, n, VariableName) \
534  BOOST_AUTO(BOOST_PP_CAT(VariableName, n), boost::fusion::at_c<BOOST_PP_CAT(index, n)>(args).GetRawPtr());
535 
536  #define NEKTAR_NEKVECTOR_UNROLL_GENERATE_VARIABLE_NAME_IN_ADDITION_SEQUENCE(z, n, VariableName) \
537  + BOOST_PP_CAT(VariableName, n)[i]
538 
539  #define NEKTAR_NEKVECTOR_UNROLL_IMPL(z, n, ClassName) \
540  template<typename IndicesType, unsigned int startIndex, unsigned int endIndex> \
541  struct ClassName<IndicesType, startIndex, endIndex, \
542  typename boost::enable_if_c \
543  < \
544  endIndex-startIndex == BOOST_PP_CAT(n, u) \
545  >::type> \
546  { \
547  BOOST_PP_REPEAT_FROM_TO(0, n, NEKTAR_NEKVECTOR_UNROLL_GENERATE_INDEX, index)\
548  \
549  template<typename AccumulatorType, typename ArgumentVectorType> \
550  static inline void Execute(AccumulatorType& accumulator, const ArgumentVectorType& args) \
551  { \
552  BOOST_AUTO(a, accumulator.GetRawPtr()); \
553  BOOST_PP_REPEAT_FROM_TO(0, n, NEKTAR_NEKVECTOR_UNROLL_GENERATE_VARIABLE, t) \
554  \
555  const unsigned int r = accumulator.GetRows(); \
556  for(unsigned int i = 0; i < r; ++i) \
557  { \
558  a[i] = t0[i] \
559  BOOST_PP_REPEAT_FROM_TO(1, n, NEKTAR_NEKVECTOR_UNROLL_GENERATE_VARIABLE_NAME_IN_ADDITION_SEQUENCE, t); \
560  } \
561  } \
562  };
563 
564  BOOST_PP_REPEAT_FROM_TO(2, NEKTAR_NEKVECTOR_MAX_UNROLL_ARGS, NEKTAR_NEKVECTOR_UNROLL_IMPL, Unroll);
565 
566  //template<typename IndicesType, unsigned int startIndex, unsigned int endIndex>
567  //struct Unroll<IndicesType, startIndex, endIndex,
568  // typename boost::enable_if_c
569  // <
570  // endIndex-startIndex == 2u
571  // >::type>
572  //{
573  // static const unsigned int index0 = boost::mpl::at_c<IndicesType, startIndex>::type::value;
574  // static const unsigned int index1 = boost::mpl::at_c<IndicesType, startIndex+1>::type::value;
575  // ;
576  // template<typename AccumulatorType, typename ArgumentVectorType>
577  // static inline void Execute(AccumulatorType& accumulator, const ArgumentVectorType& args)
578  // {
579  // BOOST_AUTO(a, accumulator.GetRawPtr());
580  // BOOST_AUTO(t0, boost::fusion::at_c<index0>(args).GetRawPtr());
581  // BOOST_AUTO(t1, boost::fusion::at_c<index1>(args).GetRawPtr());
582 
583  // const unsigned int r = accumulator.GetRows();
584  // for(unsigned int i = 0; i < r; ++i)
585  // {
586  // accumulator[i] = t0[i] + t1[i];
587  // }
588  // }
589  //};
590 
591  //template<typename IndicesType, unsigned int startIndex, unsigned int endIndex>
592  //struct Unroll<IndicesType, startIndex, endIndex,
593  // typename boost::enable_if_c
594  // <
595  // endIndex-startIndex == 3u
596  // >::type>
597  //{
598  // static const unsigned int index0 = boost::mpl::at_c<IndicesType, startIndex>::type::value;
599  // static const unsigned int index1 = boost::mpl::at_c<IndicesType, startIndex+1>::type::value;
600  // static const unsigned int index2 = boost::mpl::at_c<IndicesType, startIndex+2>::type::value;
601 
602  // template<typename AccumulatorType, typename ArgumentVectorType>
603  // static inline void Execute(AccumulatorType& accumulator, const ArgumentVectorType& args)
604  // {
605  // BOOST_AUTO(a, accumulator.GetRawPtr());
606  // BOOST_AUTO(t0, boost::fusion::at_c<index0>(args).GetRawPtr());
607  // BOOST_AUTO(t1, boost::fusion::at_c<index1>(args).GetRawPtr());
608  // BOOST_AUTO(t2, boost::fusion::at_c<index2>(args).GetRawPtr());
609 
610  // const unsigned int r = accumulator.GetRows();
611  // for(unsigned int i = 0; i < r; ++i)
612  // {
613  // accumulator[i] = t0[i] + t1[i] + t2[i];
614  // }
615  // }
616  //};
617 
618  //template<typename IndicesType, unsigned int startIndex, unsigned int endIndex>
619  //struct Unroll<IndicesType, startIndex, endIndex,
620  // typename boost::enable_if_c
621  // <
622  // endIndex-startIndex == 4u
623  // >::type>
624  //{
625  // static const unsigned int index0 = boost::mpl::at_c<IndicesType, startIndex>::type::value;
626  // static const unsigned int index1 = boost::mpl::at_c<IndicesType, startIndex+1>::type::value;
627  // static const unsigned int index2 = boost::mpl::at_c<IndicesType, startIndex+2>::type::value;
628  // static const unsigned int index3 = boost::mpl::at_c<IndicesType, startIndex+3>::type::value;
629 
630  // template<typename AccumulatorType, typename ArgumentVectorType>
631  // static inline void Execute(AccumulatorType& accumulator, const ArgumentVectorType& args)
632  // {
633  // BOOST_AUTO(a, accumulator.GetRawPtr());
634  // BOOST_AUTO(t0, boost::fusion::at_c<index0>(args).GetRawPtr());
635  // BOOST_AUTO(t1, boost::fusion::at_c<index1>(args).GetRawPtr());
636  // BOOST_AUTO(t2, boost::fusion::at_c<index2>(args).GetRawPtr());
637  // BOOST_AUTO(t3, boost::fusion::at_c<index3>(args).GetRawPtr());
638 
639  // const unsigned int r = accumulator.GetRows();
640  // for(unsigned int i = 0; i < r; ++i)
641  // {
642  // accumulator[i] = t0[i] + t1[i] + t2[i] + t3[i];
643  // }
644  // }
645  //};
646 
647  // Conditions
648  // Lhs and Rhs must result in a vector.
649  // Op must be Plus or Minus
650  // Must apply recursively.
651  template<typename LhsType, typename Op, typename RhsType, typename IndicesType, unsigned int index>
652  struct BinaryBinaryEvaluateNodeOverride<LhsType, Op, RhsType, IndicesType, index,
653  typename boost::enable_if
654  <
655  NodeCanUnroll<expt::Node<LhsType, Op, RhsType> >
656  >::type> : public boost::true_type
657  {
658  static const int endIndex = index + LhsType::TotalCount + RhsType::TotalCount;
659 
660  template<typename ResultType, typename ArgumentVectorType>
661  static inline void Evaluate(ResultType& accumulator, const ArgumentVectorType& args)
662  {
663  Unroll<IndicesType, index, endIndex>::Execute(accumulator, args);
664  }
665  };
666 }
667 #endif //NEKTAR_USE_EXPRESSION_TEMPLATES
668 
669 
670 #endif // NEKTAR_LIB_UTILITIES_NEK_VECTOR_HPP
671 
Generated on Sun Mar 15 2015 00:11:33 for Nektar++ by   doxygen 1.8.1.2