libMesh
dense_vector.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_DENSE_VECTOR_H
21 #define LIBMESH_DENSE_VECTOR_H
22 
23 // Local Includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/compare_types.h"
26 #include "libmesh/dense_vector_base.h"
27 #include "libmesh/tensor_tools.h"
28 
29 #ifdef LIBMESH_HAVE_EIGEN
30 #include "libmesh/ignore_warnings.h"
31 #include <Eigen/Core>
32 #include "libmesh/restore_warnings.h"
33 #endif
34 
35 // C++ includes
36 #include <vector>
37 
38 namespace libMesh
39 {
40 
51 template<typename T>
52 class DenseVector : public DenseVectorBase<T>
53 {
54 public:
55 
59  explicit
60  DenseVector(const unsigned int n=0);
61 
66  explicit
67  DenseVector(const unsigned int n,
68  const T & val);
69 
73  template <typename T2>
74  DenseVector (const DenseVector<T2> & other_vector);
75 
79  template <typename T2>
80  DenseVector (const std::vector<T2> & other_vector);
81 
86 
87  virtual unsigned int size() const libmesh_override
88  {
89  return cast_int<unsigned int>(_val.size());
90  }
91 
92  virtual bool empty() const libmesh_override
93  { return _val.empty(); }
94 
95  virtual void zero() libmesh_override;
96 
100  const T & operator() (const unsigned int i) const;
101 
105  T & operator() (const unsigned int i);
106 
107  virtual T el(const unsigned int i) const libmesh_override
108  { return (*this)(i); }
109 
110  virtual T & el(const unsigned int i) libmesh_override
111  { return (*this)(i); }
112 
118  template <typename T2>
119  DenseVector<T> & operator = (const DenseVector<T2> & other_vector);
120 
124  void swap(DenseVector<T> & other_vector);
125 
129  void resize (const unsigned int n);
130 
135  template <typename T2>
136  void append (const DenseVector<T2> & other_vector);
137 
141  void scale (const T factor);
142 
148  DenseVector<T> & operator*= (const T factor);
149 
157  template <typename T2, typename T3>
158  typename boostcopy::enable_if_c<
159  ScalarTraits<T2>::value, void >::type
160  add (const T2 factor,
161  const DenseVector<T3> & vec);
162 
168  template <typename T2>
169  typename CompareTypes<T, T2>::supertype dot (const DenseVector<T2> & vec) const;
170 
177  template <typename T2>
179 
183  template <typename T2>
184  bool operator== (const DenseVector<T2> & vec) const;
185 
189  template <typename T2>
190  bool operator!= (const DenseVector<T2> & vec) const;
191 
197  template <typename T2>
199 
205  template <typename T2>
207 
212  Real min () const;
213 
218  Real max () const;
219 
224  Real l1_norm () const;
225 
230  Real l2_norm () const;
231 
236  Real linfty_norm () const;
237 
242  void get_principal_subvector (unsigned int sub_n, DenseVector<T> & dest) const;
243 
251  std::vector<T> & get_values() { return _val; }
252 
256  const std::vector<T> & get_values() const { return _val; }
257 
258 private:
259 
263  std::vector<T> _val;
264 };
265 
266 
267 
268 // ------------------------------------------------------------
269 // DenseVector member functions
270 template<typename T>
271 inline
272 DenseVector<T>::DenseVector(const unsigned int n) :
273  _val (n, T(0.))
274 {
275 }
276 
277 
278 template<typename T>
279 inline
280 DenseVector<T>::DenseVector(const unsigned int n,
281  const T & val) :
282  _val (n, val)
283 {
284 }
285 
286 
287 
288 template<typename T>
289 template<typename T2>
290 inline
292  DenseVectorBase<T>()
293 {
294  const std::vector<T2> & other_vals = other_vector.get_values();
295 
296  _val.clear();
297 
298  const int N = cast_int<int>(other_vals.size());
299  _val.reserve(N);
300 
301  for (int i=0; i<N; i++)
302  _val.push_back(other_vals[i]);
303 }
304 
305 
306 
307 template<typename T>
308 template<typename T2>
309 inline
310 DenseVector<T>::DenseVector (const std::vector<T2> & other_vector) :
311  _val(other_vector)
312 {
313 }
314 
315 
316 
317 
318 
319 template<typename T>
320 template<typename T2>
321 inline
323 {
324  const std::vector<T2> & other_vals = other_vector.get_values();
325 
326  _val.clear();
327 
328  const int N = cast_int<int>(other_vals.size());
329  _val.reserve(N);
330 
331  for (int i=0; i<N; i++)
332  _val.push_back(other_vals[i]);
333 
334  return *this;
335 }
336 
337 
338 
339 template<typename T>
340 inline
342 {
343  _val.swap(other_vector._val);
344 }
345 
346 
347 
348 template<typename T>
349 inline
350 void DenseVector<T>::resize(const unsigned int n)
351 {
352  _val.resize(n);
353 
354  zero();
355 }
356 
357 
358 
359 template<typename T>
360 template<typename T2>
361 inline
362 void DenseVector<T>::append (const DenseVector<T2> & other_vector)
363 {
364  const std::vector<T2> & other_vals = other_vector.get_values();
365 
366  _val.reserve(this->size() + other_vals.size());
367  _val.insert(_val.end(), other_vals.begin(), other_vals.end());
368 }
369 
370 
371 
372 template<typename T>
373 inline
375 {
376  std::fill (_val.begin(),
377  _val.end(),
378  T(0.));
379 }
380 
381 
382 
383 template<typename T>
384 inline
385 const T & DenseVector<T>::operator () (const unsigned int i) const
386 {
387  libmesh_assert_less (i, _val.size());
388 
389  return _val[i];
390 }
391 
392 
393 
394 template<typename T>
395 inline
396 T & DenseVector<T>::operator () (const unsigned int i)
397 {
398  libmesh_assert_less (i, _val.size());
399 
400  return _val[i];
401 }
402 
403 
404 
405 template<typename T>
406 inline
407 void DenseVector<T>::scale (const T factor)
408 {
409  const int N = cast_int<int>(_val.size());
410  for (int i=0; i<N; i++)
411  _val[i] *= factor;
412 }
413 
414 
415 
416 template<typename T>
417 inline
419 {
420  this->scale(factor);
421  return *this;
422 }
423 
424 
425 
426 template<typename T>
427 template<typename T2, typename T3>
428 inline
429 typename boostcopy::enable_if_c<
430  ScalarTraits<T2>::value, void >::type
431 DenseVector<T>::add (const T2 factor,
432  const DenseVector<T3> & vec)
433 {
434  libmesh_assert_equal_to (this->size(), vec.size());
435 
436  const int N = cast_int<int>(_val.size());
437  for (int i=0; i<N; i++)
438  (*this)(i) += static_cast<T>(factor)*vec(i);
439 }
440 
441 
442 
443 template<typename T>
444 template<typename T2>
445 inline
447 {
448  if (!_val.size())
449  return 0.;
450 
451  libmesh_assert_equal_to (this->size(), vec.size());
452 
453 #ifdef LIBMESH_HAVE_EIGEN
454  // We reverse the order of the arguments to dot() here since
455  // the convention in Eigen is to take the complex conjugate of the
456  // *first* argument, while ours is to take the complex conjugate of
457  // the second.
458  return Eigen::Map<const typename Eigen::Matrix<T2, Eigen::Dynamic, 1>>(&(vec.get_values()[0]), vec.size())
459  .dot(Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(&_val[0], _val.size()));
460 #else
461  typename CompareTypes<T, T2>::supertype val = 0.;
462 
463  const int N = cast_int<int>(_val.size());
464  // The following pragma tells clang's vectorizer that it is safe to
465  // reorder floating point operations for this loop.
466 #ifdef __clang__
467 #pragma clang loop vectorize(enable)
468 #endif
469  for (int i=0; i<N; i++)
470  val += (*this)(i)*libmesh_conj(vec(i));
471 
472  return val;
473 #endif
474 }
475 
476 template<typename T>
477 template<typename T2>
478 inline
480 {
481  libmesh_assert_equal_to (this->size(), vec.size());
482 
483  typename CompareTypes<T, T2>::supertype val = 0.;
484 
485  const int N = cast_int<int>(_val.size());
486  for (int i=0; i<N; i++)
487  val += (*this)(i)*(vec(i));
488 
489  return val;
490 }
491 
492 template<typename T>
493 template<typename T2>
494 inline
496 {
497  libmesh_assert_equal_to (this->size(), vec.size());
498 
499  const int N = cast_int<int>(_val.size());
500  for (int i=0; i<N; i++)
501  if ((*this)(i) != vec(i))
502  return false;
503 
504  return true;
505 }
506 
507 
508 
509 template<typename T>
510 template<typename T2>
511 inline
513 {
514  libmesh_assert_equal_to (this->size(), vec.size());
515 
516  const int N = cast_int<int>(_val.size());
517  for (int i=0; i<N; i++)
518  if ((*this)(i) != vec(i))
519  return true;
520 
521  return false;
522 }
523 
524 
525 
526 template<typename T>
527 template<typename T2>
528 inline
530 {
531  libmesh_assert_equal_to (this->size(), vec.size());
532 
533  const int N = cast_int<int>(_val.size());
534  for (int i=0; i<N; i++)
535  (*this)(i) += vec(i);
536 
537  return *this;
538 }
539 
540 
541 
542 template<typename T>
543 template<typename T2>
544 inline
546 {
547  libmesh_assert_equal_to (this->size(), vec.size());
548 
549  const int N = cast_int<int>(_val.size());
550  for (int i=0; i<N; i++)
551  (*this)(i) -= vec(i);
552 
553  return *this;
554 }
555 
556 
557 
558 template<typename T>
559 inline
561 {
562  libmesh_assert (this->size());
563  Real my_min = libmesh_real((*this)(0));
564 
565  const int N = cast_int<int>(_val.size());
566  for (int i=1; i!=N; i++)
567  {
568  Real current = libmesh_real((*this)(i));
569  my_min = (my_min < current? my_min : current);
570  }
571  return my_min;
572 }
573 
574 
575 
576 template<typename T>
577 inline
579 {
580  libmesh_assert (this->size());
581  Real my_max = libmesh_real((*this)(0));
582 
583  const int N = cast_int<int>(_val.size());
584  for (int i=1; i!=N; i++)
585  {
586  Real current = libmesh_real((*this)(i));
587  my_max = (my_max > current? my_max : current);
588  }
589  return my_max;
590 }
591 
592 
593 
594 template<typename T>
595 inline
597 {
598  if (!_val.size())
599  return 0.;
600 
601 #ifdef LIBMESH_HAVE_EIGEN
602  return Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(&_val[0], _val.size()).template lpNorm<1>();
603 #else
604  Real my_norm = 0.;
605  const int N = cast_int<int>(_val.size());
606  for (int i=0; i!=N; i++)
607  my_norm += std::abs((*this)(i));
608 
609  return my_norm;
610 #endif
611 }
612 
613 
614 
615 template<typename T>
616 inline
618 {
619  if (!_val.size())
620  return 0.;
621 
622 #ifdef LIBMESH_HAVE_EIGEN
623  return Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(&_val[0], _val.size()).norm();
624 #else
625  Real my_norm = 0.;
626  const int N = cast_int<int>(_val.size());
627  // The following pragma tells clang's vectorizer that it is safe to
628  // reorder floating point operations for this loop.
629 #ifdef __clang__
630 #pragma clang loop vectorize(enable)
631 #endif
632  for (int i=0; i!=N; i++)
633  my_norm += TensorTools::norm_sq((*this)(i));
634 
635  return sqrt(my_norm);
636 #endif
637 }
638 
639 
640 
641 template<typename T>
642 inline
644 {
645  if (!_val.size())
646  return 0.;
647 
648 #ifdef LIBMESH_HAVE_EIGEN
649  return Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(&_val[0], _val.size()).template lpNorm<Eigen::Infinity>();
650 #else
651  Real my_norm = TensorTools::norm_sq((*this)(0));
652 
653  const int N = cast_int<int>(_val.size());
654  for (int i=1; i!=N; i++)
655  {
656  Real current = TensorTools::norm_sq((*this)(i));
657  my_norm = (my_norm > current? my_norm : current);
658  }
659  return sqrt(my_norm);
660 #endif
661 }
662 
663 
664 
665 template<typename T>
666 inline
668  DenseVector<T> & dest) const
669 {
670  libmesh_assert_less_equal ( sub_n, this->size() );
671 
672  dest.resize(sub_n);
673  const int N = cast_int<int>(sub_n);
674  for (int i=0; i<N; i++)
675  dest(i) = _val[i];
676 }
677 
678 } // namespace libMesh
679 
680 #endif // LIBMESH_DENSE_VECTOR_H
T libmesh_real(T a)
DenseVector< T > & operator*=(const T factor)
Multiplies every element in the vector by factor.
Definition: dense_vector.h:418
double abs(double a)
boostcopy::enable_if_c< ScalarTraits< T2 >::value, void >::type add(const T2 factor, const DenseVector< T3 > &vec)
Adds factor times vec to this vector.
Definition: dense_vector.h:431
T libmesh_conj(T a)
virtual T el(const unsigned int i) const libmesh_override
Definition: dense_vector.h:107
void resize(const unsigned int n)
Resize the vector.
Definition: dense_vector.h:350
CompareTypes< T, T2 >::supertype dot(const DenseVector< T2 > &vec) const
Definition: dense_vector.h:446
virtual void zero() libmesh_override
Set every element in the vector to 0.
Definition: dense_vector.h:374
DenseVector< T > & operator=(const DenseVector< T2 > &other_vector)
Assignment operator.
Definition: dense_vector.h:322
std::vector< T > & get_values()
Definition: dense_vector.h:251
Real min() const
Definition: dense_vector.h:560
Real l2_norm() const
Definition: dense_vector.h:617
const T & operator()(const unsigned int i) const
Definition: dense_vector.h:385
The libMesh namespace provides an interface to certain functionality in the library.
void scale(const T factor)
Multiplies every element in the vector by factor.
Definition: dense_vector.h:407
virtual unsigned int size() const libmesh_override
Definition: dense_vector.h:87
virtual bool empty() const libmesh_override
Definition: dense_vector.h:92
void swap(DenseVector< T > &other_vector)
STL-like swap method.
Definition: dense_vector.h:341
libmesh_assert(j)
bool operator==(const DenseVector< T2 > &vec) const
Definition: dense_vector.h:495
const std::vector< T > & get_values() const
Definition: dense_vector.h:256
Real max() const
Definition: dense_vector.h:578
Real linfty_norm() const
Definition: dense_vector.h:643
Real l1_norm() const
Definition: dense_vector.h:596
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:263
bool operator!=(const DenseVector< T2 > &vec) const
Definition: dense_vector.h:512
~DenseVector()
Destructor.
Definition: dense_vector.h:85
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual T & el(const unsigned int i) libmesh_override
Definition: dense_vector.h:110
CompareTypes< T, T2 >::supertype indefinite_dot(const DenseVector< T2 > &vec) const
Definition: dense_vector.h:479
DenseVector(const unsigned int n=0)
Constructor.
Definition: dense_vector.h:272
void append(const DenseVector< T2 > &other_vector)
Append additional entries to (resizing, but unchanging) the vector.
Definition: dense_vector.h:362
Defines an abstract dense vector base class for use in Finite Element-type computations.
Definition: dof_map.h:62
Defines a dense vector for use in Finite Element-type computations.
DenseVector< T > & operator+=(const DenseVector< T2 > &vec)
Adds vec to this vector.
Definition: dense_vector.h:529
void get_principal_subvector(unsigned int sub_n, DenseVector< T > &dest) const
Puts the principal subvector of size sub_n (i.e.
Definition: dense_vector.h:667
DenseVector< T > & operator-=(const DenseVector< T2 > &vec)
Subtracts vec from this vector.
Definition: dense_vector.h:545