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 <meshless_interpolation_function.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)
 
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
 
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 39 of file meshless_interpolation_function.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 289 of file dense_vector.h.

289  :
290  _val (n, T{})
291 {
292 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ 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 297 of file dense_vector.h.

298  :
299  _val (n, val)
300 {
301 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ DenseVector() [3/7]

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

Copy-constructor.

Definition at line 308 of file dense_vector.h.

308  :
309  DenseVectorBase<T>()
310 {
311  const std::vector<T2> & other_vals = other_vector.get_values();
312 
313  _val.clear();
314 
315  const int N = cast_int<int>(other_vals.size());
316  _val.reserve(N);
317 
318  for (int i=0; i<N; i++)
319  _val.push_back(other_vals[i]);
320 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ 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 327 of file dense_vector.h.

327  :
328  _val(other_vector)
329 {
330 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ 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 336 of file dense_vector.h.

336  :
337  _val(init_list.begin(), init_list.end())
338 {
339 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ 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 455 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().

457 {
458  libmesh_assert_equal_to (this->size(), vec.size());
459 
460  const int N = cast_int<int>(_val.size());
461  for (int i=0; i<N; i++)
462  (*this)(i) += static_cast<T>(factor)*vec(i);
463 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
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 386 of file dense_vector.h.

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

387 {
388  const std::vector<T2> & other_vals = other_vector.get_values();
389 
390  _val.reserve(this->size() + other_vals.size());
391  _val.insert(_val.end(), other_vals.begin(), other_vals.end());
392 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ 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 470 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().

471 {
472  if (!_val.size())
473  return 0.;
474 
475  libmesh_assert_equal_to (this->size(), vec.size());
476 
477 #ifdef LIBMESH_HAVE_EIGEN
478  // We reverse the order of the arguments to dot() here since
479  // the convention in Eigen is to take the complex conjugate of the
480  // *first* argument, while ours is to take the complex conjugate of
481  // the second.
482  return Eigen::Map<const typename Eigen::Matrix<T2, Eigen::Dynamic, 1>>(vec.get_values().data(), vec.size())
483  .dot(Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(_val.data(), _val.size()));
484 #else
485  typename CompareTypes<T, T2>::supertype val = 0.;
486 
487  const int N = cast_int<int>(_val.size());
488  // The following pragma tells clang's vectorizer that it is safe to
489  // reorder floating point operations for this loop.
490 #ifdef __clang__
491 #pragma clang loop vectorize(enable)
492 #endif
493  for (int i=0; i<N; i++)
494  val += (*this)(i)*libmesh_conj(vec(i));
495 
496  return val;
497 #endif
498 }
T libmesh_conj(T a)
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
CompareTypes< T, T2 >::supertype dot(const DenseVector< T2 > &vec) const
Definition: dense_vector.h:470
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 124 of file dense_vector.h.

125  { 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 127 of file dense_vector.h.

128  { 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:280

◆ 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 691 of file dense_vector.h.

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

693 {
694  libmesh_assert_less_equal ( sub_n, this->size() );
695 
696  dest.resize(sub_n);
697  const int N = cast_int<int>(sub_n);
698  for (int i=0; i<N; i++)
699  dest(i) = _val[i];
700 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
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 268 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().

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

◆ 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 273 of file dense_vector.h.

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

◆ 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 503 of file dense_vector.h.

504 {
505  libmesh_assert_equal_to (this->size(), vec.size());
506 
507  typename CompareTypes<T, T2>::supertype val = 0.;
508 
509  const int N = cast_int<int>(_val.size());
510  for (int i=0; i<N; i++)
511  val += (*this)(i)*(vec(i));
512 
513  return val;
514 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
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 620 of file dense_vector.h.

621 {
622  if (!_val.size())
623  return 0.;
624 
625 #ifdef LIBMESH_HAVE_EIGEN
626  return Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(_val.data(), _val.size()).template lpNorm<1>();
627 #else
628  Real my_norm = 0.;
629  const int N = cast_int<int>(_val.size());
630  for (int i=0; i!=N; i++)
631  my_norm += std::abs((*this)(i));
632 
633  return my_norm;
634 #endif
635 }
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:57
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
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 641 of file dense_vector.h.

642 {
643  if (!_val.size())
644  return 0.;
645 
646 #ifdef LIBMESH_HAVE_EIGEN
647  return Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(_val.data(), _val.size()).norm();
648 #else
649  Real my_norm = 0.;
650  const int N = cast_int<int>(_val.size());
651  // The following pragma tells clang's vectorizer that it is safe to
652  // reorder floating point operations for this loop.
653 #ifdef __clang__
654 #pragma clang loop vectorize(enable)
655 #endif
656  for (int i=0; i!=N; i++)
657  my_norm += TensorTools::norm_sq((*this)(i));
658 
659  return sqrt(my_norm);
660 #endif
661 }
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:53
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
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 667 of file dense_vector.h.

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

668 {
669  if (!_val.size())
670  return 0.;
671 
672 #ifdef LIBMESH_HAVE_EIGEN
673  return Eigen::Map<const typename Eigen::Matrix<T, Eigen::Dynamic, 1>>(_val.data(), _val.size()).template lpNorm<Eigen::Infinity>();
674 #else
675  Real my_norm = TensorTools::norm_sq((*this)(0));
676 
677  const int N = cast_int<int>(_val.size());
678  for (int i=1; i!=N; i++)
679  {
680  Real current = TensorTools::norm_sq((*this)(i));
681  my_norm = (my_norm > current? my_norm : current);
682  }
683  return sqrt(my_norm);
684 #endif
685 }
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:53
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
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 602 of file dense_vector.h.

603 {
604  libmesh_assert (this->size());
605  Real my_max = libmesh_real((*this)(0));
606 
607  const int N = cast_int<int>(_val.size());
608  for (int i=1; i!=N; i++)
609  {
610  Real current = libmesh_real((*this)(i));
611  my_max = (my_max > current? my_max : current);
612  }
613  return my_max;
614 }
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:280
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 584 of file dense_vector.h.

585 {
586  libmesh_assert (this->size());
587  Real my_min = libmesh_real((*this)(0));
588 
589  const int N = cast_int<int>(_val.size());
590  for (int i=1; i!=N; i++)
591  {
592  Real current = libmesh_real((*this)(i));
593  my_min = (my_min < current? my_min : current);
594  }
595  return my_min;
596 }
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:280
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 536 of file dense_vector.h.

537 {
538  libmesh_assert_equal_to (this->size(), vec.size());
539 
540  const int N = cast_int<int>(_val.size());
541  for (int i=0; i<N; i++)
542  if ((*this)(i) != vec(i))
543  return true;
544 
545  return false;
546 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
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 409 of file dense_vector.h.

410 {
411  libmesh_assert_less (i, _val.size());
412 
413  return _val[i];
414 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ 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 420 of file dense_vector.h.

421 {
422  libmesh_assert_less (i, _val.size());
423 
424  return _val[i];
425 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ 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 442 of file dense_vector.h.

443 {
444  this->scale(factor);
445  return *this;
446 }
void scale(const T factor)
Multiplies every element in the vector by factor.
Definition: dense_vector.h:431

◆ 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 553 of file dense_vector.h.

554 {
555  libmesh_assert_equal_to (this->size(), vec.size());
556 
557  const int N = cast_int<int>(_val.size());
558  for (int i=0; i<N; i++)
559  (*this)(i) += vec(i);
560 
561  return *this;
562 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
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 569 of file dense_vector.h.

570 {
571  libmesh_assert_equal_to (this->size(), vec.size());
572 
573  const int N = cast_int<int>(_val.size());
574  for (int i=0; i<N; i++)
575  (*this)(i) -= vec(i);
576 
577  return *this;
578 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
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 346 of file dense_vector.h.

347 {
348  const std::vector<T2> & other_vals = other_vector.get_values();
349 
350  _val.clear();
351 
352  const int N = cast_int<int>(other_vals.size());
353  _val.reserve(N);
354 
355  for (int i=0; i<N; i++)
356  _val.push_back(other_vals[i]);
357 
358  return *this;
359 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ 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 519 of file dense_vector.h.

520 {
521  libmesh_assert_equal_to (this->size(), vec.size());
522 
523  const int N = cast_int<int>(_val.size());
524  for (int i=0; i<N; i++)
525  if ((*this)(i) != vec(i))
526  return false;
527 
528  return true;
529 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280
virtual unsigned int size() const override final
Definition: dense_vector.h:104

◆ print()

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

Pretty-print the vector to stdout.

Definition at line 51 of file dense_vector_base.C.

References libMesh::make_range().

52 {
53  for (auto i : make_range(this->size()))
54  os << std::setw(8)
55  << this->el(i)
56  << std::endl;
57 }
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:134
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:134
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 374 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_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(), libMesh::DGFEMContext::neighbor_side_fe_reinit(), SkewFunc::operator()(), BdyFunction::operator()(), libMesh::MeshlessInterpolationFunction::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(), 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().

375 {
376  _val.resize(n);
377 
378  zero();
379 }
virtual void zero() override final
Set every element in the vector to 0.
Definition: dense_vector.h:398
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ scale()

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

Multiplies every element in the vector by factor.

Definition at line 431 of file dense_vector.h.

Referenced by SolidSystem::element_time_derivative().

432 {
433  const int N = cast_int<int>(_val.size());
434  for (int i=0; i<N; i++)
435  _val[i] *= factor;
436 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ 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()(), SlitFunc::operator()(), PeriodicQuadFunction::operator()(), SolutionFunction< dim >::operator()(), libMesh::ConstFunction< Output >::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:280

◆ swap()

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

STL-like swap method.

Definition at line 365 of file dense_vector.h.

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

366 {
367  _val.swap(other_vector._val);
368 }
std::vector< T > _val
The actual data values, stored as a 1D array.
Definition: dense_vector.h:280

◆ 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 actual data values, stored as a 1D array.

Definition at line 280 of file dense_vector.h.

Referenced by libMesh::DenseVector< Output >::empty(), libMesh::DenseVector< Output >::get_values(), and libMesh::DenseVector< Output >::size().


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