Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
NekVector.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File: NekVector.cpp
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:
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
38 
39 
40 namespace Nektar
41 {
42 
43  template<typename DataType>
45  m_size(0),
46  m_data(),
47  m_wrapperType(eCopy)
48  {
49  }
50 
51  template<typename DataType>
52  NekVector<DataType>::NekVector(unsigned int size) :
53  m_size(size),
54  m_data(size),
55  m_wrapperType(eCopy)
56  {
57  }
58 
59  template<typename DataType>
60  NekVector<DataType>::NekVector(unsigned int size, typename boost::call_traits<DataType>::const_reference a) :
61  m_size(size),
62  m_data(size),
63  m_wrapperType(eCopy)
64  {
65  std::fill_n(m_data.get(), m_size, a);
66  }
67 
68  template<typename DataType>
69  NekVector<DataType>::NekVector(const std::string& vectorValues) :
70  m_size(0),
71  m_data(),
72  m_wrapperType(eCopy)
73  {
74  try
75  {
76  std::vector<DataType> values = FromString<DataType>(vectorValues);
77  m_size = values.size();
79  std::copy(values.begin(), values.end(), m_data.begin());
80 
81  ASSERTL0(m_size > 0, "Error converting string values to vector");
82  }
83  catch(std::runtime_error& e)
84  {
85  NEKERROR(ErrorUtil::efatal, e.what());
86  }
87  }
88 
89  template<typename DataType>
90  NekVector<DataType>::NekVector(typename boost::call_traits<DataType>::const_reference x,
91  typename boost::call_traits<DataType>::const_reference y,
92  typename boost::call_traits<DataType>::const_reference z) :
93  m_size(3),
94  m_data(m_size),
95  m_wrapperType(eCopy)
96  {
97  m_data[0] = x;
98  m_data[1] = y;
99  m_data[2] = z;
100  }
101 
102  template<typename DataType>
104  m_size(rhs.GetDimension()),
105  m_data(rhs.m_data),
106  m_wrapperType(rhs.m_wrapperType)
107  {
108  if( m_wrapperType == eCopy )
109  {
111  std::copy(rhs.begin(), rhs.end(), m_data.get());
112  }
113  }
114 
115  template<typename DataType>
116  NekVector<DataType>::NekVector(unsigned int size, const DataType* const ptr) :
117  m_size(size),
118  m_data(size, ptr),
119  m_wrapperType(eCopy)
120  {
121  }
122 
123  template<typename DataType>
125  m_size(ptr.num_elements()),
126  m_data(ptr),
127  m_wrapperType(h)
128  {
129  if( h == eCopy )
130  {
132  CopyArray(ptr, m_data);
133  }
134  }
135 
136  template<typename DataType>
138  m_size(size),
139  m_data(ptr),
140  m_wrapperType(h)
141  {
142  if( h == eCopy )
143  {
144  ASSERTL0(size <= ptr.num_elements(), "Attempting to populate a vector of size " +
145  boost::lexical_cast<std::string>(size) + " but the incoming array only has " +
146  boost::lexical_cast<std::string>(ptr.num_elements()) + " elements.");
147 
149  std::copy(ptr.begin(), ptr.begin()+size, m_data.begin());
150  }
151  }
152 
153  template<typename DataType>
155  m_size(ptr.num_elements()),
156  m_data(ptr, eVECTOR_WRAPPER),
157  m_wrapperType(h)
158  {
159  if( h == eCopy )
160  {
162  CopyArray(ptr, m_data);
163  }
164  }
165 
166  template<typename DataType>
168  m_size(size),
169  m_data(ptr, eVECTOR_WRAPPER),
170  m_wrapperType(h)
171  {
172  if( h == eCopy )
173  {
174  ASSERTL0(size <= ptr.num_elements(), "Attempting to populate a vector of size " +
175  boost::lexical_cast<std::string>(size) + " but the incoming array only has " +
176  boost::lexical_cast<std::string>(ptr.num_elements()) + " elements.");
177 
179  std::copy(ptr.begin(), ptr.begin()+size, m_data.begin());
180  }
181  }
182 
183  template<typename DataType>
185 
186  template<typename DataType>
188  {
189  if( m_wrapperType == eCopy )
190  {
191  // If the current vector is a copy, then regardless of the rhs type
192  // we just copy over the values, resizing if needed.
193  if( GetDimension() != rhs.GetDimension() )
194  {
195  m_size = rhs.GetDimension();
196  m_data = Array<OneD, DataType>(m_size);
197  }
198  }
199  else if( m_wrapperType == eWrapper )
200  {
201  // If the current vector is wrapped, then just copy over the top,
202  // but the sizes of the two vectors must be the same.
203  ASSERTL0(GetDimension() == rhs.GetDimension(), "Wrapped NekVectors must have the same dimension in operator=");
204  }
205 
206  std::copy(rhs.begin(), rhs.end(), m_data.get());
207  return *this;
208  }
209 
210 
211  template<typename DataType>
213  {
214  return m_size;
215  }
216 
217  template<typename DataType>
218  unsigned int NekVector<DataType>::GetRows() const
219  {
220  return m_size;
221  }
222 
223  template<typename DataType>
225  {
226  return this->GetData().get();
227  }
228 
229  template<typename DataType>
230  Array<OneD, DataType>& NekVector<DataType>::GetPtr() { return this->GetData(); }
231 
232  template<typename DataType>
233  const DataType* NekVector<DataType>::GetRawPtr() const
234  {
235  return m_data.get();
236  }
237 
238  template<typename DataType>
240 
241  template<typename DataType>
243 
244  template<typename DataType>
245  typename NekVector<DataType>::iterator NekVector<DataType>::end() { return GetRawPtr() + this->GetDimension(); }
246 
247  template<typename DataType>
248  typename NekVector<DataType>::const_iterator NekVector<DataType>::begin() const { return GetRawPtr(); }
249 
250  template<typename DataType>
251  typename NekVector<DataType>::const_iterator NekVector<DataType>::end() const { return GetRawPtr() + GetDimension(); }
252 
253  template<typename DataType>
254  typename boost::call_traits<DataType>::reference NekVector<DataType>::operator()(unsigned int i)
255  {
256  ASSERTL1((i >= 0) && (i < this->GetDimension()), "Invalid access to m_data via parenthesis operator");
257  return this->GetData()[i];
258  }
259 
260  template<typename DataType>
261  typename boost::call_traits<DataType>::reference NekVector<DataType>::operator[](unsigned int i)
262  {
263  return this->GetData()[i];
264  }
265 
266  template<typename DataType>
267  typename boost::call_traits<DataType>::reference NekVector<DataType>::x()
268  {
269  ASSERTL1(this->GetDimension() >= 1, "Invalid use of NekVector::x");
270  return (*this)(0);
271  }
272 
273  template<typename DataType>
274  typename boost::call_traits<DataType>::reference NekVector<DataType>::y()
275  {
276  ASSERTL1(this->GetDimension() >= 2, "Invalid use of NekVector::y");
277  return (*this)(1);
278  }
279 
280  template<typename DataType>
281  typename boost::call_traits<DataType>::reference NekVector<DataType>::z()
282  {
283  ASSERTL1(this->GetDimension() >= 3, "Invalid use of NekVector::z");
284  return (*this)(2);
285  }
286 
287  template<typename DataType>
288  void NekVector<DataType>::SetX(typename boost::call_traits<DataType>::const_reference val)
289  {
290  ASSERTL1(this->GetDimension() >= 1, "Invalid use of NekVector::SetX");
291  this->GetData()[0] = val;
292  }
293 
294  template<typename DataType>
295  void NekVector<DataType>::SetY(typename boost::call_traits<DataType>::const_reference val)
296  {
297  ASSERTL1(this->GetDimension() >= 2, "Invalid use of NekVector::SetX");
298  this->GetData()[1] = val;
299  }
300 
301  template<typename DataType>
302  void NekVector<DataType>::SetZ(typename boost::call_traits<DataType>::const_reference val)
303  {
304  ASSERTL1(this->GetDimension() >= 3, "Invalid use of NekVector::SetX");
305  this->GetData()[2] = val;
306  }
307 
308  template<typename DataType>
310  {
311  AddEqual(*this, rhs);
312  return *this;
313  }
314 
315  template<typename DataType>
317  {
318  SubtractEqual(*this, rhs);
319  return *this;
320  }
321 
322  template<typename DataType>
323  NekVector<DataType>& NekVector<DataType>::operator*=(typename boost::call_traits<DataType>::const_reference rhs)
324  {
325  MultiplyEqual(*this, rhs);
326  return *this;
327  }
328 
329  template<typename DataType>
330  NekVector<DataType>& NekVector<DataType>::operator/=(typename boost::call_traits<DataType>::const_reference rhs)
331  {
332  DivideEqual(*this, rhs);
333  return *this;
334  }
335 
336  template<typename DataType>
338 
339  template<typename DataType>
340  typename boost::call_traits<DataType>::const_reference NekVector<DataType>::operator()(unsigned int i) const
341  {
342  ASSERTL1(( i >= 0) && (i < GetDimension()), "Invalid access to m_data via parenthesis operator");
343  return m_data[i];
344  }
345 
346  template<typename DataType>
347  typename boost::call_traits<DataType>::const_reference NekVector<DataType>::operator[](unsigned int i) const
348  {
349  return m_data[i];
350  }
351 
352  template<typename DataType>
353  typename boost::call_traits<DataType>::const_reference NekVector<DataType>::x() const
354  {
355  ASSERTL1( GetDimension() >= 1, "Invalid use of NekVector::x");
356  return (*this)(0);
357  }
358 
359  template<typename DataType>
360  typename boost::call_traits<DataType>::const_reference NekVector<DataType>::y() const
361  {
362  ASSERTL1( GetDimension() >= 2, "Invalid use of NekVector::y");
363  return (*this)(1);
364  }
365 
366  template<typename DataType>
367  typename boost::call_traits<DataType>::const_reference NekVector<DataType>::z() const
368  {
369  ASSERTL1( GetDimension() >= 3, "Invalid use of NekVector::z");
370  return (*this)(2);
371  }
372 
373  #ifndef NEKTAR_USE_EXPRESSION_TEMPLATES
374  template<typename DataType>
376  #endif
377 
378  template<typename DataType>
379  DataType NekVector<DataType>::Magnitude() const { return Nektar::Magnitude(*this); }
380 
381  template<typename DataType>
382  DataType NekVector<DataType>::Dot(const NekVector<DataType>& rhs) const { return Nektar::Dot(*this, rhs); }
383 
384  template<typename DataType>
386  {
387  return Nektar::Cross(*this, rhs);
388  }
389 
390  template<typename DataType>
391  std::string NekVector<DataType>::AsString() const { return Nektar::AsString(*this); }
392 
393  // Norms
394  template<typename DataType>
395  DataType NekVector<DataType>::L1Norm() const { return Nektar::L1Norm(*this); }
396 
397  template<typename DataType>
398  DataType NekVector<DataType>::L2Norm() const { return Nektar::L2Norm(*this); }
399 
400  template<typename DataType>
401  DataType NekVector<DataType>::InfinityNorm() const { return Nektar::InfinityNorm(*this); }
402 
403  template<typename DataType>
404  PointerWrapper NekVector<DataType>::GetWrapperType() const { return m_wrapperType; }
405 
406  template<typename DataType>
408 
409  template<typename DataType>
410  void NekVector<DataType>::SetSize(unsigned int s) { m_size = s; }
411 
412  template<typename DataType>
414 
415  template<typename DataType>
416  void NekVector<DataType>::SetData(const Array<OneD, DataType>& newData) { m_data = newData; }
417 
418  template<typename DataType>
419  void NekVector<DataType>::Resize(unsigned int newSize)
420  {
421  if(m_data.num_elements() < newSize )
422  {
423  m_data = Array<OneD, DataType>(newSize);
424  }
425  m_size = newSize;
426  }
427 
429 
430  template<typename DataType>
431  void Add(NekVector<DataType>& result,
432  const NekVector<DataType>& lhs,
433  const NekVector<DataType>& rhs)
434  {
435  DataType* r_buf = result.GetRawPtr();
436  const DataType* lhs_buf = lhs.GetRawPtr();
437  const DataType* rhs_buf = rhs.GetRawPtr();
438  const unsigned int ldim = lhs.GetDimension();
439  for(int i = 0; i < ldim; ++i)
440  {
441  r_buf[i] = lhs_buf[i] + rhs_buf[i];
442  }
443  }
444 
445  template<typename DataType>
447  const NekVector<DataType>& lhs,
448  const NekVector<DataType>& rhs)
449  {
450  DataType* r_buf = result.GetRawPtr();
451  const DataType* lhs_buf = lhs.GetRawPtr();
452  const DataType* rhs_buf = rhs.GetRawPtr();
453  const unsigned int ldim = lhs.GetDimension();
454  for(int i = 0; i < ldim; ++i)
455  {
456  r_buf[i] = -lhs_buf[i] + rhs_buf[i];
457  }
458  }
459 
460  template LIB_UTILITIES_EXPORT void Add(NekVector<NekDouble>& result,
461  const NekVector<NekDouble>& lhs,
462  const NekVector<NekDouble>& rhs);
463  template LIB_UTILITIES_EXPORT void AddNegatedLhs(NekVector<NekDouble>& result,
464  const NekVector<NekDouble>& lhs,
465  const NekVector<NekDouble>& rhs);
466 
467  template<typename DataType>
469  const NekVector<DataType>& rhs)
470  {
471  DataType* r_buf = result.GetRawPtr();
472  const DataType* rhs_buf = rhs.GetRawPtr();
473  const unsigned int rdim = rhs.GetDimension();
474  for(int i = 0; i < rdim; ++i)
475  {
476  r_buf[i] += rhs_buf[i];
477  }
478  }
479 
480  template<typename DataType>
482  const NekVector<DataType>& rhs)
483  {
484  DataType* r_buf = result.GetRawPtr();
485  const DataType* rhs_buf = rhs.GetRawPtr();
486  const unsigned int rdim = rhs.GetDimension();
487  for(int i = 0; i < rdim; ++i)
488  {
489  r_buf[i] = -r_buf[i] + rhs_buf[i];
490  }
491  }
492 
493  template LIB_UTILITIES_EXPORT
494  void AddEqual(NekVector<NekDouble>& result,
495  const NekVector<NekDouble>& rhs);
496  template LIB_UTILITIES_EXPORT
497  void AddEqualNegatedLhs(NekVector<NekDouble>& result,
498  const NekVector<NekDouble>& rhs);
499 
500  template<typename LhsDataType,
501  typename RhsDataType>
503  const NekVector<RhsDataType>& rhs)
504  {
505  NekVector<LhsDataType> result(lhs.GetDimension());
506  Add(result, lhs, rhs);
507  return result;
508  }
509 
510  template LIB_UTILITIES_EXPORT
511  NekVector<NekDouble> Add(const NekVector<NekDouble>& lhs,
512  const NekVector<NekDouble>& rhs);
513 
514  template<typename ResultDataType, typename InputDataType>
517  const NekVector<InputDataType>& rhs)
518  {
519  ResultDataType* r_buf = result.GetRawPtr();
520  typename boost::add_const<InputDataType>::type* lhs_buf = lhs.GetRawPtr();
521  typename boost::add_const<InputDataType>::type* rhs_buf = rhs.GetRawPtr();
522  const unsigned int ldim = lhs.GetDimension();
523  for(int i = 0; i < ldim; ++i)
524  {
525  r_buf[i] = lhs_buf[i] - rhs_buf[i];
526  }
527  }
528 
529  template<typename ResultDataType, typename InputDataType>
532  const NekVector<InputDataType>& rhs)
533  {
534  ResultDataType* r_buf = result.GetRawPtr();
535  typename boost::add_const<InputDataType>::type* lhs_buf = lhs.GetRawPtr();
536  typename boost::add_const<InputDataType>::type* rhs_buf = rhs.GetRawPtr();
537  const unsigned int ldim = lhs.GetDimension();
538  for(int i = 0; i < ldim; ++i)
539  {
540  r_buf[i] = -lhs_buf[i] - rhs_buf[i];
541  }
542  }
543 
544  template LIB_UTILITIES_EXPORT
545  void Subtract(NekVector<NekDouble>& result,
546  const NekVector<NekDouble>& lhs,
547  const NekVector<NekDouble>& rhs);
548 
549  template LIB_UTILITIES_EXPORT
550  void SubtractNegatedLhs(NekVector<NekDouble>& result,
551  const NekVector<NekDouble>& lhs,
552  const NekVector<NekDouble>& rhs);
553 
554  template<typename ResultDataType, typename InputDataType>
556  const NekVector<InputDataType>& rhs)
557  {
558  ResultDataType* r_buf = result.GetRawPtr();
559  typename boost::add_const<InputDataType>::type* rhs_buf = rhs.GetRawPtr();
560  const unsigned int rdim = rhs.GetDimension();
561  for(int i = 0; i < rdim; ++i)
562  {
563  r_buf[i] -= rhs_buf[i];
564  }
565  }
566 
567  template<typename ResultDataType, typename InputDataType>
569  const NekVector<InputDataType>& rhs)
570  {
571  ResultDataType* r_buf = result.GetRawPtr();
572  typename boost::add_const<InputDataType>::type* rhs_buf = rhs.GetRawPtr();
573  const unsigned int rdim = rhs.GetDimension();
574  for(int i = 0; i < rdim; ++i)
575  {
576  r_buf[i] = -r_buf[i] - rhs_buf[i];
577  }
578  }
579 
580  template LIB_UTILITIES_EXPORT
581  void SubtractEqual(NekVector<NekDouble>& result,
582  const NekVector<NekDouble>& rhs);
583 
584  template LIB_UTILITIES_EXPORT
585  void SubtractEqualNegatedLhs(NekVector<NekDouble>& result,
586  const NekVector<NekDouble>& rhs);
587 
588  template<typename DataType>
589  NekVector<DataType>
591  const NekVector<DataType>& rhs)
592  {
593  NekVector<DataType> result(lhs.GetDimension());
594  Subtract(result, lhs, rhs);
595  return result;
596  }
597 
598  template LIB_UTILITIES_EXPORT
599  NekVector<NekDouble>
600  Subtract(const NekVector<NekDouble>& lhs,
601  const NekVector<NekDouble>& rhs);
602 
603  template<typename ResultDataType, typename InputDataType>
606  const NekDouble& rhs)
607  {
608  ResultDataType* r_buf = result.GetRawPtr();
609  typename boost::add_const<InputDataType>::type* lhs_buf = lhs.GetRawPtr();
610 
611  const unsigned int ldim = lhs.GetDimension();
612  for(int i = 0; i < ldim; ++i)
613  {
614  r_buf[i] = lhs_buf[i] / rhs;
615  }
616  }
617 
618  template LIB_UTILITIES_EXPORT
619  void Divide(NekVector<NekDouble>& result,
620  const NekVector<NekDouble>& lhs,
621  const NekDouble& rhs);
622 
623  template<typename ResultDataType>
625  const NekDouble& rhs)
626  {
627  ResultDataType* r_buf = result.GetRawPtr();
628 
629  const unsigned int resdim = result.GetDimension();
630  for(int i = 0; i < resdim; ++i)
631  {
632  r_buf[i] /= rhs;
633  }
634  }
635 
636  template LIB_UTILITIES_EXPORT
637  void DivideEqual(NekVector<NekDouble>& result,
638  const NekDouble& rhs);
639 
640  template<typename DataType>
641  NekVector<DataType>
643  const NekDouble& rhs)
644  {
645  NekVector<DataType> result(lhs.GetDimension());
646  Divide(result, lhs, rhs);
647  return result;
648  }
649 
650  template LIB_UTILITIES_EXPORT
651  NekVector<NekDouble>
652  Divide(const NekVector<NekDouble>& lhs,
653  const NekDouble& rhs);
654 
655 
656  template<typename ResultDataType, typename InputDataType>
659  const NekVector<InputDataType>& rhs)
660  {
661  ResultDataType* result_buf = result.GetRawPtr();
662  const InputDataType* rhs_buf = rhs.GetRawPtr();
663  const InputDataType* lhs_buf = lhs.GetRawPtr();
664  const unsigned int resdim = result.GetDimension();
665  for(int i = 0; i < resdim; ++i)
666  {
667  result_buf[i] = lhs_buf[i] * rhs_buf[i];
668  }
669  }
670 
671  template LIB_UTILITIES_EXPORT
672  void Multiply(NekVector<NekDouble>& result, const NekVector<NekDouble>& lhs, const NekVector<NekDouble>& rhs);
673 
674  template<typename ResultDataType, typename InputDataType>
676  const NekVector<InputDataType>& rhs)
677  {
678  ResultDataType* result_buf = result.GetRawPtr();
679  const InputDataType* rhs_buf = rhs.GetRawPtr();
680  const unsigned int resdim = result.GetDimension();
681  for(int i = 0; i < resdim; ++i)
682  {
683  result_buf[i] *= rhs_buf[i];
684  }
685  }
686 
687  template LIB_UTILITIES_EXPORT
688  void MultiplyEqual(NekVector<NekDouble>& result, const NekVector<NekDouble>& rhs);
689 
690  template<typename DataType, typename InputDataType>
691  NekVector<DataType>
693  const NekVector<InputDataType>& rhs)
694  {
695  NekVector<DataType> result(lhs.GetDimension());
696  Multiply(result, lhs, rhs);
697  return result;
698  }
699 
700  template LIB_UTILITIES_EXPORT
701  NekVector<NekDouble>
702  Multiply(const NekVector<NekDouble>& lhs, const NekVector<NekDouble>& rhs);
703 
704 
705  template<typename ResultDataType, typename InputDataType>
708  const NekDouble& rhs)
709  {
710  ResultDataType* r_buf = result.GetRawPtr();
711  const InputDataType* lhs_buf = lhs.GetRawPtr();
712 
713  const unsigned int ldim = lhs.GetDimension();
714  for(int i = 0; i < ldim; ++i)
715  {
716  r_buf[i] = lhs_buf[i] * rhs;
717  }
718  }
719 
720  template LIB_UTILITIES_EXPORT
721  void Multiply(NekVector<NekDouble>& result,
722  const NekVector<NekDouble>& lhs,
723  const NekDouble& rhs);
724 
725  template<typename ResultDataType>
727  const NekDouble& rhs)
728  {
729  ResultDataType* r_buf = result.GetRawPtr();
730  const unsigned int rdim = result.GetDimension();
731  for(unsigned int i = 0; i < rdim; ++i)
732  {
733  r_buf[i] *= rhs;
734  }
735  }
736  template LIB_UTILITIES_EXPORT
737  void MultiplyEqual(NekVector<NekDouble>& result,
738  const NekDouble& rhs);
739 
740  template<typename DataType>
741  NekVector<DataType>
743  const NekDouble& rhs)
744  {
745  NekVector<DataType> result(lhs.GetDimension());
746  Multiply(result, lhs, rhs);
747  return result;
748  }
749 
750  template LIB_UTILITIES_EXPORT
751  NekVector<NekDouble>
752  Multiply(const NekVector<NekDouble>& lhs,
753  const NekDouble& rhs);
754 
755  template<typename ResultDataType, typename InputDataType>
757  const NekDouble& lhs,
758  const NekVector<InputDataType>& rhs)
759  {
760  Multiply(result, rhs, lhs);
761  }
762 
763  template<typename ResultDataType, typename InputDataType>
765  const NekDouble& lhs,
766  const NekVector<InputDataType>& rhs)
767  {
768  ResultDataType* r_buf = result.GetRawPtr();
769  const InputDataType* rhs_buf = rhs.GetRawPtr();
770  NekDouble inverse = 1.0/lhs;
771 
772  const unsigned int rdim = rhs.GetDimension();
773  for(int i = 0; i < rdim; ++i)
774  {
775  r_buf[i] = inverse * rhs_buf[i];
776  }
777  }
778 
779  template LIB_UTILITIES_EXPORT
780  void MultiplyInvertedLhs(NekVector<NekDouble>& result,
781  const NekDouble& lhs,
782  const NekVector<NekDouble>& rhs);
783 
784  template LIB_UTILITIES_EXPORT
785  void Multiply(NekVector<NekDouble>& result,
786  const NekDouble& lhs,
787  const NekVector<NekDouble>& rhs);
788 
789  template<typename DataType>
791  const NekVector<DataType>& rhs)
792  {
793  return Multiply(rhs, lhs);
794  }
795 
796  template LIB_UTILITIES_EXPORT
797  NekVector<NekDouble> Multiply(const NekDouble& lhs,
798  const NekVector<NekDouble>& rhs);
799 
800 
801  template<typename DataType>
802  std::ostream& operator<<(std::ostream& os, const NekVector<DataType>& rhs)
803  {
804  os << rhs.AsString();
805  return os;
806  }
807 
808  template LIB_UTILITIES_EXPORT
809  std::ostream& operator<<(std::ostream& os, const NekVector<NekDouble>& rhs);
810 
811  template<typename DataType>
813  const NekPoint<DataType>& dest)
814  {
815  NekVector<DataType> result(3, 0.0);
816  for(unsigned int i = 0; i < 3; ++i)
817  {
818  result[i] = dest[i]-source[i];
819  }
820  return result;
821  }
822 
823 
824  template LIB_UTILITIES_EXPORT
825  NekVector<NekDouble> createVectorFromPoints(const NekPoint<NekDouble>& source,
826  const NekPoint<NekDouble>& dest);
827 
828  template<typename DataType>
830  const DataType& t)
831  {
832  NekPoint<DataType> result;
833  for(unsigned int i = 0; i < 3; ++i)
834  {
835  result[i] = lhs[i]*t;
836  }
837 
838  return result;
839  }
840 
841  template LIB_UTILITIES_EXPORT
842  NekPoint<NekDouble> findPointAlongVector(const NekVector<NekDouble>& lhs,
843  const NekDouble& t);
844 
845  template<typename DataType>
847  const NekVector<DataType>& rhs)
848  {
849  if( lhs.GetDimension() != rhs.GetDimension() )
850  {
851  return false;
852  }
853 
854  return std::equal(lhs.begin(), lhs.end(), rhs.begin());
855  }
856 
857  template LIB_UTILITIES_EXPORT
858  bool operator==(const NekVector<NekDouble>& lhs,
859  const NekVector<NekDouble>& rhs);
860 
861  template<typename DataType>
863  const NekVector<DataType>& rhs)
864  {
865  return !(lhs == rhs);
866  }
867 
868  template LIB_UTILITIES_EXPORT
869  bool operator!=(const NekVector<NekDouble>& lhs,
870  const NekVector<NekDouble>& rhs);
871 
872  template<typename DataType>
873  std::vector<DataType> FromString(const std::string& str)
874  {
875  std::vector<DataType> result;
876 
877  try
878  {
879  typedef boost::tokenizer<boost::char_separator<char> > tokenizer;
880  boost::char_separator<char> sep("(<,>) ");
881  tokenizer tokens(str, sep);
882  for( tokenizer::iterator strIter = tokens.begin(); strIter != tokens.end(); ++strIter)
883  {
884  result.push_back(boost::lexical_cast<DataType>(*strIter));
885  }
886  }
887  catch(boost::bad_lexical_cast&)
888  {
889  }
890 
891  return result;
892  }
893 
894  template LIB_UTILITIES_EXPORT
895  std::vector<NekDouble> FromString(const std::string& str);
896 
897  template<typename DataType>
898  DataType L1Norm(const NekVector<DataType>& v)
899  {
900  typedef NekVector<DataType> VectorType;
901 
902  DataType result(0);
903  for(typename VectorType::const_iterator iter = v.begin(); iter != v.end(); ++iter)
904  {
905  result += fabs(*iter);
906  }
907 
908  return result;
909  }
910 
911  template LIB_UTILITIES_EXPORT
912  NekDouble L1Norm(const NekVector<NekDouble>& v);
913 
914  template<typename DataType>
915  DataType L2Norm(const NekVector<DataType>& v)
916  {
917  typedef NekVector<DataType> VectorType;
918 
919  DataType result(0);
920  for(typename VectorType::const_iterator iter = v.begin(); iter != v.end(); ++iter)
921  {
922  DataType v = fabs(*iter);
923  result += v*v;
924  }
925  return sqrt(result);
926  }
927 
928  template LIB_UTILITIES_EXPORT
929  NekDouble L2Norm(const NekVector<NekDouble>& v);
930 
931  template<typename DataType>
933  {
934  DataType result = fabs(v[0]);
935  const unsigned int vdim = v.GetDimension();
936  for(unsigned int i = 0; i < vdim; ++i)
937  {
938  result = std::max(fabs(v[i]), result);
939  }
940  return result;
941  }
942 
943  template LIB_UTILITIES_EXPORT
944  NekDouble InfinityNorm(const NekVector<NekDouble>& v);
945 
946  template<typename DataType>
948  {
949  NekVector<DataType> temp(v);
950  const unsigned int tdim = temp.GetDimension();
951  for(unsigned int i = 0; i < tdim; ++i)
952  {
953  temp(i) = -temp(i);
954  }
955  return temp;
956  }
957 
958  template LIB_UTILITIES_EXPORT
959  NekVector<NekDouble> Negate(const NekVector<NekDouble>& v);
960 
961  template<typename DataType>
963  {
964  DataType* data = v.GetRawPtr();
965  const unsigned int vdim = v.GetDimension();
966  for(unsigned int i = 0; i < vdim; ++i)
967  {
968  data[i] = -data[i];
969  }
970  }
971 
972  template LIB_UTILITIES_EXPORT
973  void NegateInPlace(NekVector<NekDouble>& v);
974 
975  template<typename DataType>
976  DataType Magnitude(const NekVector<DataType>& v)
977  {
978  DataType result = DataType(0);
979 
980  const unsigned int vdim = v.GetDimension();
981  for(unsigned int i = 0; i < vdim; ++i)
982  {
983  result += v[i]*v[i];
984  }
985  return sqrt(result);
986  }
987 
988  template LIB_UTILITIES_EXPORT
989  NekDouble Magnitude(const NekVector<NekDouble>& v) ;
990 
991  template<typename DataType>
992  DataType Dot(const NekVector<DataType>& lhs,
993  const NekVector<DataType>& rhs)
994  {
995  ASSERTL1( lhs.GetDimension() == rhs.GetDimension(), "Dot, dimension of the two operands must be identical.");
996 
997  DataType result = DataType(0);
998  const unsigned int ldim = lhs.GetDimension();
999  for(unsigned int i = 0; i < ldim; ++i)
1000  {
1001  result += lhs[i]*rhs[i];
1002  }
1003 
1004  return result;
1005  }
1006 
1007  template LIB_UTILITIES_EXPORT
1008  NekDouble Dot(const NekVector<NekDouble>& lhs,
1009  const NekVector<NekDouble>& rhs) ;
1010 
1011  template<typename DataType>
1013  {
1014  DataType m = v.Magnitude();
1015  if( m > DataType(0) )
1016  {
1017  v /= m;
1018  }
1019  }
1020 
1021  template LIB_UTILITIES_EXPORT
1022  void Normalize(NekVector<NekDouble>& v);
1023 
1024  void NegateInPlace(NekDouble& v) { v = -v; }
1025  void InvertInPlace(NekDouble& v) { v = 1.0/v; }
1026 
1027  template<typename DataType>
1029  const NekVector<DataType>& rhs)
1030  {
1031  ASSERTL1(lhs.GetDimension() == 3 && rhs.GetDimension() == 3, "Cross is only valid for 3D vectors.");
1032 
1033  DataType first = lhs.y()*rhs.z() - lhs.z()*rhs.y();
1034  DataType second = lhs.z()*rhs.x() - lhs.x()*rhs.z();
1035  DataType third = lhs.x()*rhs.y() - lhs.y()*rhs.x();
1036 
1037  NekVector<DataType> result(first, second, third);
1038  return result;
1039  }
1040 
1041  template LIB_UTILITIES_EXPORT
1042  NekVector<NekDouble> Cross(const NekVector<NekDouble>& lhs, const NekVector<NekDouble>& rhs);
1043 
1044  template<typename DataType>
1045  std::string AsString(const NekVector<DataType>& v)
1046  {
1047  unsigned int d = v.GetRows();
1048  std::string result = "(";
1049  for(unsigned int i = 0; i < d; ++i)
1050  {
1051  result += boost::lexical_cast<std::string>(v[i]);
1052  if( i < v.GetDimension()-1 )
1053  {
1054  result += ", ";
1055  }
1056  }
1057  result += ")";
1058  return result;
1059  }
1060 
1061  template LIB_UTILITIES_EXPORT
1062  std::string AsString(const NekVector<NekDouble>& v);
1063 }
1064