libMesh
Public Member Functions | Private Attributes | List of all members
libMesh::DenseVector< T > Class Template Reference

Defines a dense vector for use in Finite Element-type computations. More...

#include <dof_map.h>

Inheritance diagram for libMesh::DenseVector< T >:
[legend]

Public Member Functions

 DenseVector (const unsigned int n=0)
 Constructor. More...
 
 DenseVector (const unsigned int n, const T &val)
 Constructor. More...
 
template<typename T2 >
 DenseVector (const DenseVector< T2 > &other_vector)
 Copy-constructor. More...
 
template<typename T2 >
 DenseVector (const std::vector< T2 > &other_vector)
 Copy-constructor, from a std::vector. More...
 
template<typename T2 >
 DenseVector (std::initializer_list< T2 > init_list)
 Initializer list constructor. More...
 
 DenseVector (DenseVector &&)=default
 The 5 special functions can be defaulted for this class, as it does not manage any memory itself. More...
 
 DenseVector (const DenseVector &)=default
 
DenseVectoroperator= (const DenseVector &)=default
 
DenseVectoroperator= (DenseVector &&)=default
 
virtual ~DenseVector ()=default
 
virtual unsigned int size () const override final
 
virtual bool empty () const override final
 
virtual void zero () override final
 Set every element in the vector to 0. More...
 
const T & operator() (const unsigned int i) const
 
T & operator() (const unsigned int i)
 
const T & operator[] (const unsigned int i) const
 
T & operator[] (const unsigned int i)
 
virtual T el (const unsigned int i) const override final
 
virtual T & el (const unsigned int i) override final
 
template<typename T2 >
DenseVector< T > & operator= (const DenseVector< T2 > &other_vector)
 Assignment operator. More...
 
void swap (DenseVector< T > &other_vector)
 STL-like swap method. More...
 
void resize (const unsigned int n)
 Resize the vector. More...
 
template<typename T2 >
void append (const DenseVector< T2 > &other_vector)
 Append additional entries to (resizing, but unchanging) the vector. More...
 
void scale (const T factor)
 Multiplies every element in the vector by factor. More...
 
DenseVector< T > & operator*= (const T factor)
 Multiplies every element in the vector by factor. More...
 
template<typename T2 , typename T3 >
boostcopy::enable_if_c< ScalarTraits< T2 >::value, void >::type add (const T2 factor, const DenseVector< T3 > &vec)
 Adds factor times vec to this vector. More...
 
template<typename T2 >
CompareTypes< T, T2 >::supertype dot (const DenseVector< T2 > &vec) const
 
template<typename T2 >
CompareTypes< T, T2 >::supertype indefinite_dot (const DenseVector< T2 > &vec) const
 
template<typename T2 >
bool operator== (const DenseVector< T2 > &vec) const
 
template<typename T2 >
bool operator!= (const DenseVector< T2 > &vec) const
 
template<typename T2 >
DenseVector< T > & operator+= (const DenseVector< T2 > &vec)
 Adds vec to this vector. More...
 
template<typename T2 >
DenseVector< T > & operator-= (const DenseVector< T2 > &vec)
 Subtracts vec from this vector. More...
 
Real min () const
 
Real max () const
 
Real l1_norm () const
 
Real l2_norm () const
 
Real linfty_norm () const
 
void get_principal_subvector (unsigned int sub_n, DenseVector< T > &dest) const
 Puts the principal subvector of size sub_n (i.e. More...
 
std::vector< T > & get_values ()
 
const std::vector< T > & get_values () const
 
std::vector< T >::const_iterator begin () const
 
std::vector< T >::iterator begin ()
 
std::vector< T >::const_iterator end () const
 
std::vector< T >::iterator end ()
 
void print (std::ostream &os) const
 Pretty-print the vector to stdout. More...
 
void print_scientific (std::ostream &os, unsigned precision=8) const
 Prints the entries of the vector with additional decimal places in scientific notation. More...
 

Private Attributes

std::vector< T > _val
 The actual data values, stored as a 1D array. More...
 

Detailed Description

template<typename T>
class libMesh::DenseVector< T >

Defines a dense vector for use in Finite Element-type computations.

This class is to basically compliment the DenseMatrix class. It has additional capabilities over the std::vector that make it useful for finite elements, particularly for systems of equations. All overridden virtual functions are documented in dense_vector_base.h.

Author
Benjamin S. Kirk
Date
2003

Definition at line 74 of file dof_map.h.

Constructor & Destructor Documentation

◆ DenseVector() [1/7]

template<typename T >
libMesh::DenseVector< T >::DenseVector ( const unsigned int  n = 0)
inlineexplicit

Constructor.

Creates a dense vector of dimension n.

Definition at line 311 of file dense_vector.h.

311  :
312  _val (n, T{})
313 {
314 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ DenseVector() [2/7]

template<typename T>
libMesh::DenseVector< T >::DenseVector ( const unsigned int  n,
const T &  val 
)
inlineexplicit

Constructor.

Creates a dense vector of dimension n where all entries have value val.

Definition at line 319 of file dense_vector.h.

320  :
321  _val (n, val)
322 {
323 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ DenseVector() [3/7]

template<typename T >
template<typename T2 >
libMesh::DenseVector< T >::DenseVector ( const DenseVector< T2 > &  other_vector)
inline

Copy-constructor.

Definition at line 330 of file dense_vector.h.

330  :
331  DenseVectorBase<T>()
332 {
333  const std::vector<T2> & other_vals = other_vector.get_values();
334 
335  _val.clear();
336 
337  const int N = cast_int<int>(other_vals.size());
338  _val.reserve(N);
339 
340  for (int i=0; i<N; i++)
341  _val.push_back(other_vals[i]);
342 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ DenseVector() [4/7]

template<typename T >
template<typename T2 >
libMesh::DenseVector< T >::DenseVector ( const std::vector< T2 > &  other_vector)
inline

Copy-constructor, from a std::vector.

Definition at line 349 of file dense_vector.h.

349  :
350  _val(other_vector)
351 {
352 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ DenseVector() [5/7]

template<typename T >
template<typename T2 >
libMesh::DenseVector< T >::DenseVector ( std::initializer_list< T2 >  init_list)
inline

Initializer list constructor.

Definition at line 358 of file dense_vector.h.

358  :
359  _val(init_list.begin(), init_list.end())
360 {
361 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ DenseVector() [6/7]

template<typename T>
libMesh::DenseVector< T >::DenseVector ( DenseVector< T > &&  )
default

The 5 special functions can be defaulted for this class, as it does not manage any memory itself.

◆ DenseVector() [7/7]

template<typename T>
libMesh::DenseVector< T >::DenseVector ( const DenseVector< T > &  )
default

◆ ~DenseVector()

template<typename T>
virtual libMesh::DenseVector< T >::~DenseVector ( )
virtualdefault

Member Function Documentation

◆ add()

template<typename T >
template<typename T2 , typename T3 >
boostcopy::enable_if_c< ScalarTraits< T2 >::value, void >::type libMesh::DenseVector< T >::add ( const T2  factor,
const DenseVector< T3 > &  vec 
)
inline

Adds factor times vec to this vector.

This should only work if T += T2 * T3 is valid C++ and if T2 is scalar. Return type is void

Returns
A reference to *this.

Definition at line 477 of file dense_vector.h.

Referenced by libMesh::EulerSolver::_general_residual(), libMesh::NewmarkSolver::_general_residual(), libMesh::TransientRBEvaluation::rb_solve(), libMesh::RBEvaluation::rb_solve(), libMesh::TransientRBEvaluation::rb_solve_again(), and libMesh::DenseMatrix< Real >::vector_mult_add().

479 {
480  libmesh_assert_equal_to (this->size(), vec.size());
481 
482  const int N = cast_int<int>(_val.size());
483  for (int i=0; i<N; i++)
484  (*this)(i) += static_cast<T>(factor)*vec(i);
485 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ append()

template<typename T >
template<typename T2 >
void libMesh::DenseVector< T >::append ( const DenseVector< T2 > &  other_vector)
inline

Append additional entries to (resizing, but unchanging) the vector.

Definition at line 408 of file dense_vector.h.

Referenced by libMesh::FEGenericBase< FEOutputType< T >::type >::coarsened_dof_values().

409 {
410  const std::vector<T2> & other_vals = other_vector.get_values();
411 
412  _val.reserve(this->size() + other_vals.size());
413  _val.insert(_val.end(), other_vals.begin(), other_vals.end());
414 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ begin() [1/2]

template<typename T>
std::vector<T>::const_iterator libMesh::DenseVector< T >::begin ( ) const
inline
Returns
An iterator pointing to the beginning of the vector

Definition at line 288 of file dense_vector.h.

288 { return _val.begin(); }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ begin() [2/2]

template<typename T>
std::vector<T>::iterator libMesh::DenseVector< T >::begin ( )
inline

Definition at line 289 of file dense_vector.h.

289 { return _val.begin(); }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ dot()

template<typename T >
template<typename T2 >
CompareTypes< T, T2 >::supertype libMesh::DenseVector< T >::dot ( const DenseVector< T2 > &  vec) const
inline
Returns
The dot product of *this with vec.

In the complex-valued case, uses the complex conjugate of vec.

Definition at line 492 of file dense_vector.h.

Referenced by libMesh::TransientRBConstruction::assemble_affine_expansion(), libMesh::TransientRBEvaluation::compute_residual_dual_norm(), libMesh::RBEIMEvaluation::get_eim_error_indicator(), libMesh::TransientRBEvaluation::rb_solve(), libMesh::RBEvaluation::rb_solve(), libMesh::TransientRBConstruction::set_error_temporal_data(), and libMesh::TransientRBConstruction::update_RB_initial_condition_all_N().

493 {
494  if (!_val.size())
495  return 0.;
496 
497  libmesh_assert_equal_to (this->size(), vec.size());
498 
499 #ifdef LIBMESH_HAVE_EIGEN
500  // We reverse the order of the arguments to dot() here since
501  // the convention in Eigen is to take the complex conjugate of the
502  // *first* argument, while ours is to take the complex conjugate of
503  // the second.
504  return Eigen::Map<const typename Eigen::Matrix<T2, Eigen::Dynamic, 1>>(vec.get_values().data(), vec.size())
505  .dot(Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(_val.data(), _val.size()));
506 #else
507  typename CompareTypes<T, T2>::supertype val = 0.;
508 
509  const int N = cast_int<int>(_val.size());
510  // The following pragma tells clang's vectorizer that it is safe to
511  // reorder floating point operations for this loop.
512 #ifdef __clang__
513 #pragma clang loop vectorize(enable)
514 #endif
515  for (int i=0; i<N; i++)
516  val += (*this)(i)*libmesh_conj(vec(i));
517 
518  return val;
519 #endif
520 }
T libmesh_conj(T a)
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
CompareTypes< T, T2 >::supertype dot(const DenseVector< T2 > &vec) const
Definition: dense_vector.h:492
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ el() [1/2]

template<typename T>
virtual T libMesh::DenseVector< T >::el ( const unsigned int  i) const
inlinefinaloverridevirtual
Returns
The (i) element of the vector.

Implements libMesh::DenseVectorBase< T >.

Definition at line 134 of file dense_vector.h.

135  { return (*this)(i); }

◆ el() [2/2]

template<typename T>
virtual T& libMesh::DenseVector< T >::el ( const unsigned int  i)
inlinefinaloverridevirtual
Returns
The (i) element of the vector as a writable reference.

Implements libMesh::DenseVectorBase< T >.

Definition at line 137 of file dense_vector.h.

138  { return (*this)(i); }

◆ empty()

template<typename T>
virtual bool libMesh::DenseVector< T >::empty ( ) const
inlinefinaloverridevirtual
Returns
true iff size() is 0.

Reimplemented from libMesh::DenseVectorBase< T >.

Definition at line 109 of file dense_vector.h.

Referenced by libMesh::NumericVector< Number >::add_vector(), libMesh::NumericVector< Number >::insert(), and MeshFunctionTest::test_subdomain_id_sets().

110  { return _val.empty(); }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ end() [1/2]

template<typename T>
std::vector<T>::const_iterator libMesh::DenseVector< T >::end ( ) const
inline
Returns
An iterator pointing to the end of the vector

Definition at line 294 of file dense_vector.h.

294 { return _val.end(); }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ end() [2/2]

template<typename T>
std::vector<T>::iterator libMesh::DenseVector< T >::end ( )
inline

Definition at line 295 of file dense_vector.h.

295 { return _val.end(); }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ get_principal_subvector()

template<typename T>
void libMesh::DenseVector< T >::get_principal_subvector ( unsigned int  sub_n,
DenseVector< T > &  dest 
) const
inline

Puts the principal subvector of size sub_n (i.e.

first sub_n entries) into dest.

Definition at line 713 of file dense_vector.h.

Referenced by libMesh::RBEIMEvaluation::get_eim_error_indicator(), and libMesh::TransientRBConstruction::update_RB_initial_condition_all_N().

715 {
716  libmesh_assert_less_equal ( sub_n, this->size() );
717 
718  dest.resize(sub_n);
719  const int N = cast_int<int>(sub_n);
720  for (int i=0; i<N; i++)
721  dest(i) = _val[i];
722 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ get_values() [1/2]

template<typename T>
std::vector<T>& libMesh::DenseVector< T >::get_values ( )
inline
Returns
A reference to the underlying data storage vector.

This should be used with caution (i.e. one should not change the size of the vector, etc.) but is useful for interoperating with low level BLAS routines which expect a simple array.

Definition at line 278 of file dense_vector.h.

Referenced by libMesh::DenseMatrix< Real >::_evd_lapack(), libMesh::DenseMatrix< Real >::_lu_back_substitute_lapack(), libMesh::DenseMatrix< Real >::_matvec_blas(), libMesh::DenseMatrix< Real >::_svd_solve_lapack(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubProjector::construct_projection(), libMesh::GradientMeshFunction::operator()(), libMesh::FEMContext::pre_fe_reinit(), and DenseMatrixTest::testEVD_helper().

278 { return _val; }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ get_values() [2/2]

template<typename T>
const std::vector<T>& libMesh::DenseVector< T >::get_values ( ) const
inline
Returns
A constant reference to the underlying data storage vector.

Definition at line 283 of file dense_vector.h.

283 { return _val; }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ indefinite_dot()

template<typename T >
template<typename T2 >
CompareTypes< T, T2 >::supertype libMesh::DenseVector< T >::indefinite_dot ( const DenseVector< T2 > &  vec) const
inline
Returns
The dot product of *this with vec.

In the complex-valued case, does not use the complex conjugate of vec.

Definition at line 525 of file dense_vector.h.

526 {
527  libmesh_assert_equal_to (this->size(), vec.size());
528 
529  typename CompareTypes<T, T2>::supertype val = 0.;
530 
531  const int N = cast_int<int>(_val.size());
532  for (int i=0; i<N; i++)
533  val += (*this)(i)*(vec(i));
534 
535  return val;
536 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ l1_norm()

template<typename T >
Real libMesh::DenseVector< T >::l1_norm ( ) const
inline
Returns
The \(l_1\)-norm of the vector, i.e. the sum of the absolute values of the entries.

Definition at line 642 of file dense_vector.h.

643 {
644  if (!_val.size())
645  return 0.;
646 
647 #ifdef LIBMESH_HAVE_EIGEN
648  return Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(_val.data(), _val.size()).template lpNorm<1>();
649 #else
650  Real my_norm = 0.;
651  const int N = cast_int<int>(_val.size());
652  for (int i=0; i!=N; i++)
653  my_norm += std::abs((*this)(i));
654 
655  return my_norm;
656 #endif
657 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ l2_norm()

template<typename T >
Real libMesh::DenseVector< T >::l2_norm ( ) const
inline
Returns
The \(l_2\)-norm of the vector, i.e. the square root of the sum of the squares of the entries.

Definition at line 663 of file dense_vector.h.

664 {
665  if (!_val.size())
666  return 0.;
667 
668 #ifdef LIBMESH_HAVE_EIGEN
669  return Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(_val.data(), _val.size()).norm();
670 #else
671  Real my_norm = 0.;
672  const int N = cast_int<int>(_val.size());
673  // The following pragma tells clang's vectorizer that it is safe to
674  // reorder floating point operations for this loop.
675 #ifdef __clang__
676 #pragma clang loop vectorize(enable)
677 #endif
678  for (int i=0; i!=N; i++)
679  my_norm += TensorTools::norm_sq((*this)(i));
680 
681  return sqrt(my_norm);
682 #endif
683 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
auto norm(const T &a) -> decltype(std::abs(a))
Definition: tensor_tools.h:74
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
auto norm_sq(const T &a) -> decltype(std::norm(a))
Definition: tensor_tools.h:104

◆ linfty_norm()

template<typename T >
Real libMesh::DenseVector< T >::linfty_norm ( ) const
inline
Returns
The \(l_\infty\)-norm of the vector, i.e. the maximum absolute value of the entries.

Definition at line 689 of file dense_vector.h.

Referenced by libMesh::RBEIMEvaluation::get_eim_error_indicator().

690 {
691  if (!_val.size())
692  return 0.;
693 
694 #ifdef LIBMESH_HAVE_EIGEN
695  return Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(_val.data(), _val.size()).template lpNorm<Eigen::Infinity>();
696 #else
697  Real my_norm = TensorTools::norm_sq((*this)(0));
698 
699  const int N = cast_int<int>(_val.size());
700  for (int i=1; i!=N; i++)
701  {
702  Real current = TensorTools::norm_sq((*this)(i));
703  my_norm = (my_norm > current? my_norm : current);
704  }
705  return sqrt(my_norm);
706 #endif
707 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
auto norm_sq(const T &a) -> decltype(std::norm(a))
Definition: tensor_tools.h:104

◆ max()

template<typename T >
Real libMesh::DenseVector< T >::max ( ) const
inline
Returns
The maximum entry of the vector, or the maximum real part in the case of complex numbers.

Definition at line 624 of file dense_vector.h.

625 {
626  libmesh_assert (this->size());
627  Real my_max = libmesh_real((*this)(0));
628 
629  const int N = cast_int<int>(_val.size());
630  for (int i=1; i!=N; i++)
631  {
632  Real current = libmesh_real((*this)(i));
633  my_max = (my_max > current? my_max : current);
634  }
635  return my_max;
636 }
T libmesh_real(T a)
libmesh_assert(ctx)
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ min()

template<typename T >
Real libMesh::DenseVector< T >::min ( ) const
inline
Returns
The minimum entry of the vector, or the minimum real part in the case of complex numbers.

Definition at line 606 of file dense_vector.h.

607 {
608  libmesh_assert (this->size());
609  Real my_min = libmesh_real((*this)(0));
610 
611  const int N = cast_int<int>(_val.size());
612  for (int i=1; i!=N; i++)
613  {
614  Real current = libmesh_real((*this)(i));
615  my_min = (my_min < current? my_min : current);
616  }
617  return my_min;
618 }
T libmesh_real(T a)
libmesh_assert(ctx)
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ operator!=()

template<typename T >
template<typename T2 >
bool libMesh::DenseVector< T >::operator!= ( const DenseVector< T2 > &  vec) const
inline
Returns
true if vec is not exactly equal to this vector, false otherwise.

Definition at line 558 of file dense_vector.h.

559 {
560  libmesh_assert_equal_to (this->size(), vec.size());
561 
562  const int N = cast_int<int>(_val.size());
563  for (int i=0; i<N; i++)
564  if ((*this)(i) != vec(i))
565  return true;
566 
567  return false;
568 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ operator()() [1/2]

template<typename T >
const T & libMesh::DenseVector< T >::operator() ( const unsigned int  i) const
inline
Returns
Entry i of the vector as a const reference.

Definition at line 431 of file dense_vector.h.

432 {
433  libmesh_assert_less (i, _val.size());
434 
435  return _val[i];
436 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ operator()() [2/2]

template<typename T >
T & libMesh::DenseVector< T >::operator() ( const unsigned int  i)
inline
Returns
Entry i of the vector as a writable reference.

Definition at line 442 of file dense_vector.h.

443 {
444  libmesh_assert_less (i, _val.size());
445 
446  return _val[i];
447 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ operator*=()

template<typename T>
DenseVector< T > & libMesh::DenseVector< T >::operator*= ( const T  factor)
inline

Multiplies every element in the vector by factor.

Returns
A reference to *this.

Definition at line 464 of file dense_vector.h.

465 {
466  this->scale(factor);
467  return *this;
468 }
void scale(const T factor)
Multiplies every element in the vector by factor.
Definition: dense_vector.h:453

◆ operator+=()

template<typename T >
template<typename T2 >
DenseVector< T > & libMesh::DenseVector< T >::operator+= ( const DenseVector< T2 > &  vec)
inline

Adds vec to this vector.

Returns
A reference to *this.

Definition at line 575 of file dense_vector.h.

576 {
577  libmesh_assert_equal_to (this->size(), vec.size());
578 
579  const int N = cast_int<int>(_val.size());
580  for (int i=0; i<N; i++)
581  (*this)(i) += vec(i);
582 
583  return *this;
584 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ operator-=()

template<typename T >
template<typename T2 >
DenseVector< T > & libMesh::DenseVector< T >::operator-= ( const DenseVector< T2 > &  vec)
inline

Subtracts vec from this vector.

Returns
A reference to *this.

Definition at line 591 of file dense_vector.h.

592 {
593  libmesh_assert_equal_to (this->size(), vec.size());
594 
595  const int N = cast_int<int>(_val.size());
596  for (int i=0; i<N; i++)
597  (*this)(i) -= vec(i);
598 
599  return *this;
600 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ operator=() [1/3]

template<typename T>
DenseVector& libMesh::DenseVector< T >::operator= ( const DenseVector< T > &  )
default

◆ operator=() [2/3]

template<typename T>
DenseVector& libMesh::DenseVector< T >::operator= ( DenseVector< T > &&  )
default

◆ operator=() [3/3]

template<typename T >
template<typename T2 >
DenseVector< T > & libMesh::DenseVector< T >::operator= ( const DenseVector< T2 > &  other_vector)
inline

Assignment operator.

Returns
A reference to *this.

Definition at line 368 of file dense_vector.h.

369 {
370  const std::vector<T2> & other_vals = other_vector.get_values();
371 
372  _val.clear();
373 
374  const int N = cast_int<int>(other_vals.size());
375  _val.reserve(N);
376 
377  for (int i=0; i<N; i++)
378  _val.push_back(other_vals[i]);
379 
380  return *this;
381 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ operator==()

template<typename T >
template<typename T2 >
bool libMesh::DenseVector< T >::operator== ( const DenseVector< T2 > &  vec) const
inline
Returns
true if vec is exactly equal to this vector, false otherwise.

Definition at line 541 of file dense_vector.h.

542 {
543  libmesh_assert_equal_to (this->size(), vec.size());
544 
545  const int N = cast_int<int>(_val.size());
546  for (int i=0; i<N; i++)
547  if ((*this)(i) != vec(i))
548  return false;
549 
550  return true;
551 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ operator[]() [1/2]

template<typename T>
const T& libMesh::DenseVector< T >::operator[] ( const unsigned int  i) const
inline
Returns
Entry i of the vector as a const reference.

Definition at line 127 of file dense_vector.h.

127 { return (*this)(i); }

◆ operator[]() [2/2]

template<typename T>
T& libMesh::DenseVector< T >::operator[] ( const unsigned int  i)
inline
Returns
Entry i of the vector as a writable reference.

Definition at line 132 of file dense_vector.h.

132 { return (*this)(i); }

◆ print()

template<typename T >
void libMesh::DenseVectorBase< T >::print ( std::ostream &  os) const
inherited

Pretty-print the vector to stdout.

Definition at line 110 of file dense_vector_base.h.

References libMesh::make_range().

111 {
112  for (auto i : make_range(this->size()))
113  os << std::setw(8)
114  << this->el(i)
115  << std::endl;
116 }
virtual unsigned int size() const =0
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:140
virtual T el(const unsigned int i) const =0

◆ print_scientific()

template<typename T >
void libMesh::DenseVectorBase< T >::print_scientific ( std::ostream &  os,
unsigned  precision = 8 
) const
inherited

Prints the entries of the vector with additional decimal places in scientific notation.

Definition at line 31 of file dense_vector_base.C.

References libMesh::make_range().

32 {
33  // save the initial format flags
34  std::ios_base::fmtflags os_flags = os.flags();
35 
36  // Print the vector entries.
37  for (auto i : make_range(this->size()))
38  os << std::setw(10)
39  << std::scientific
40  << std::setprecision(precision)
41  << this->el(i)
42  << std::endl;
43 
44  // reset the original format flags
45  os.flags(os_flags);
46 }
virtual unsigned int size() const =0
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:140
virtual T el(const unsigned int i) const =0

◆ resize()

template<typename T >
void libMesh::DenseVector< T >::resize ( const unsigned int  n)
inline

Resize the vector.

Sets all elements to 0.

Definition at line 396 of file dense_vector.h.

Referenced by libMesh::DenseMatrix< Real >::_cholesky_back_substitute(), libMesh::DenseMatrix< Real >::_evd_lapack(), libMesh::DenseMatrix< Real >::_lu_back_substitute(), libMesh::DenseMatrix< Real >::_svd_lapack(), libMesh::DenseMatrix< Real >::_svd_solve_lapack(), libMesh::HPCoarsenTest::add_projection(), libMesh::TransientRBConstruction::allocate_data_structures(), alternative_fe_assembly(), assemble(), LinearElasticity::assemble(), assemble_1D(), AssembleOptimization::assemble_A_and_F(), assemble_biharmonic(), assemble_cd(), assemble_divgrad(), assemble_elasticity(), assemble_ellipticdg(), assemble_func(), assemble_graddiv(), assemble_helmholtz(), assemble_laplace(), assemble_matrix_and_rhs(), assemble_poisson(), assemble_shell(), assemble_stokes(), assemble_wave(), Biharmonic::JR::bounds(), libMesh::DofMap::build_constraint_matrix_and_vector(), libMesh::TransientRBEvaluation::cache_online_residual_terms(), libMesh::FEGenericBase< FEOutputType< T >::type >::coarsened_dof_values(), compute_enriched_soln(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_periodic_constraints(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_proj_constraints(), compute_residual(), SolidSystem::element_time_derivative(), fe_assembly(), form_functionA(), form_functionB(), libMesh::DGFEMContext::neighbor_side_fe_reinit(), SkewFunc::operator()(), BdyFunction::operator()(), libMesh::AutoAreaFunction::operator()(), libMesh::MeshfreeInterpolationFunction::operator()(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::CompositeFEMFunction< Output >::operator()(), libMesh::CompositeFunction< Output >::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::MeshFunction::operator()(), libMesh::MeshTools::Generation::Private::GaussLobattoRedistributionFunction::operator()(), libMesh::BoundaryProjectSolution::operator()(), periodic_bc_test_poisson(), libMesh::FEMContext::pre_fe_reinit(), libMesh::TransientRBEvaluation::rb_solve(), libMesh::RBEvaluation::rb_solve(), libMesh::TransientRBEvaluation::rb_solve_again(), libMesh::HDGProblem::residual(), LaplaceYoung::residual(), LargeDeformationElasticity::residual(), Biharmonic::JR::residual_and_jacobian(), LinearElasticityWithContact::residual_and_jacobian(), libMesh::DenseVector< Output >::resize(), libMesh::TransientRBEvaluation::resize_data_structures(), libMesh::HPCoarsenTest::select_refinement(), libMesh::DenseMatrix< Real >::vector_mult(), libMesh::DenseMatrix< Real >::vector_mult_add(), and libMesh::DenseMatrix< Real >::vector_mult_transpose().

397 {
398  _val.resize(n);
399 
400  zero();
401 }
virtual void zero() override final
Set every element in the vector to 0.
Definition: dense_vector.h:420
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ scale()

template<typename T>
void libMesh::DenseVector< T >::scale ( const T  factor)
inline

Multiplies every element in the vector by factor.

Definition at line 453 of file dense_vector.h.

Referenced by SolidSystem::element_time_derivative().

454 {
455  const int N = cast_int<int>(_val.size());
456  for (int i=0; i<N; i++)
457  _val[i] *= factor;
458 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ size()

template<typename T>
virtual unsigned int libMesh::DenseVector< T >::size ( ) const
inlinefinaloverridevirtual
Returns
The size of the vector.

Implements libMesh::DenseVectorBase< T >.

Definition at line 104 of file dense_vector.h.

Referenced by libMesh::EulerSolver::_general_residual(), libMesh::Euler2Solver::_general_residual(), libMesh::NewmarkSolver::_general_residual(), libMesh::MeshFunction::_gradient_on_elem(), libMesh::DenseMatrix< Real >::_lu_back_substitute(), libMesh::DenseMatrix< Real >::_matvec_blas(), libMesh::DenseMatrix< Real >::_svd_lapack(), libMesh::DenseMatrix< Real >::_svd_solve_lapack(), libMesh::HPCoarsenTest::add_projection(), libMesh::RBDataSerialization::add_rb_eim_evaluation_data_to_builder(), libMesh::NumericVector< Number >::add_vector(), libMesh::FEMSystem::assembly(), libMesh::DofMap::constrain_element_residual(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubProjector::construct_projection(), libMesh::RBEIMEvaluation::decrement_vector(), libMesh::MeshFunction::discontinuous_value(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::RBEIMEvaluation::get_eim_error_indicator(), libMesh::MeshFunction::hessian(), libMesh::DofMap::heterogeneously_constrain_element_jacobian_and_residual(), libMesh::DofMap::heterogeneously_constrain_element_residual(), libMesh::VectorSetAction< Val >::insert(), libMesh::NumericVector< Number >::insert(), libMesh::LaplacianErrorEstimator::internal_side_integration(), libMesh::RBConstruction::load_rb_solution(), libMesh::RBEIMEvaluation::node_decrement_vector(), ExampleOneFunction::operator()(), PeriodicQuadFunction::operator()(), SlitFunc::operator()(), SolutionFunction< dim >::operator()(), libMesh::ConstFunction< Output >::operator()(), libMesh::FDMGradient< GradType >::operator()(), libMesh::WrappedFunction< Output >::operator()(), libMesh::ParsedFEMFunction< T >::operator()(), libMesh::CompositeFEMFunction< Output >::operator()(), libMesh::CompositeFunction< Output >::operator()(), libMesh::MeshFunction::operator()(), TripleFunction::operator()(), libMesh::DenseMatrix< Real >::outer_product(), periodic_bc_test_poisson(), libMesh::RBEIMEvaluation::rb_eim_solve(), libMesh::TransientRBEvaluation::rb_solve_again(), libMesh::RBEIMEvaluation::side_decrement_vector(), DenseMatrixTest::testComplexSVD(), DenseMatrixTest::testEVD_helper(), DenseMatrixTest::testOuterProduct(), DenseMatrixTest::testSVD(), MetaPhysicL::RawType< libMesh::DenseVector< T > >::value(), libMesh::DenseMatrix< Real >::vector_mult(), libMesh::DenseMatrix< Real >::vector_mult_add(), and libMesh::DenseMatrix< Real >::vector_mult_transpose().

105  {
106  return cast_int<unsigned int>(_val.size());
107  }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ swap()

template<typename T>
void libMesh::DenseVector< T >::swap ( DenseVector< T > &  other_vector)
inline

STL-like swap method.

Definition at line 387 of file dense_vector.h.

Referenced by libMesh::EulerSolver::_general_residual(), libMesh::Euler2Solver::_general_residual(), and libMesh::DenseVector< Output >::swap().

388 {
389  _val.swap(other_vector._val);
390 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:302

◆ zero()

template<typename T >
void libMesh::DenseVector< T >::zero ( )
inlinefinaloverridevirtual

Member Data Documentation

◆ _val

template<typename T>
std::vector<T> libMesh::DenseVector< T >::_val
private

The documentation for this class was generated from the following files: