libMesh::DenseVector< T > Class Template Reference

#include <analytic_function.h>

Inheritance diagram for libMesh::DenseVector< T >:

Public Member Functions

 DenseVector (const unsigned int n=0)
 
template<typename T2 >
 DenseVector (const DenseVector< T2 > &other_vector)
 
template<typename T2 >
 DenseVector (const std::vector< T2 > &other_vector)
 
 ~DenseVector ()
 
virtual unsigned int size () const
 
virtual void zero ()
 
operator() (const unsigned int i) const
 
T & operator() (const unsigned int i)
 
virtual T el (const unsigned int i) const
 
virtual T & el (const unsigned int i)
 
template<typename T2 >
DenseVector< T > & operator= (const DenseVector< T2 > &other_vector)
 
void swap (DenseVector< T > &other_vector)
 
void resize (const unsigned int n)
 
void scale (const T factor)
 
DenseVector< T > & operator*= (const T factor)
 
template<typename T2 , typename T3 >
boostcopy::enable_if_c
< ScalarTraits< T2 >::value,
void >::type 
add (const T2 factor, const DenseVector< T3 > &vec)
 
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)
 
template<typename T2 >
DenseVector< T > & operator-= (const DenseVector< T2 > &vec)
 
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
 
std::vector< T > & get_values ()
 
const std::vector< T > & get_values () const
 
void print (std::ostream &os) const
 
void print_scientific (std::ostream &os) const
 

Private Attributes

std::vector< T > _val
 

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 DenseMatix class. It has additional capabilities over the std::vector that make it useful for finite elements, particulary for systems of equations.

Author
Benjamin S. Kirk, 2003

Definition at line 36 of file analytic_function.h.

Constructor & Destructor Documentation

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

Constructor. Creates a dense vector of dimension n.

Definition at line 254 of file dense_vector.h.

254  :
255  _val (n, T(0.))
256 {
257 }
template<typename T >
template<typename T2 >
libMesh::DenseVector< T >::DenseVector ( const DenseVector< T2 > &  other_vector)
inline

Copy-constructor.

Definition at line 264 of file dense_vector.h.

References libMesh::DenseVector< T >::_val.

264  :
265  DenseVectorBase<T>()
266 {
267  const std::vector<T2> &other_vals = other_vector.get_values();
268 
269  _val.clear();
270  _val.reserve(other_vals.size());
271 
272  for (unsigned int i=0; i<other_vals.size(); i++)
273  _val.push_back(other_vals[i]);
274 }
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 281 of file dense_vector.h.

281  :
282  _val(other_vector)
283 {
284 }
template<typename T>
libMesh::DenseVector< T >::~DenseVector ( )
inline

Destructor. Does nothing.

Definition at line 76 of file dense_vector.h.

76 {}

Member Function Documentation

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

Definition at line 388 of file dense_vector.h.

Referenced by libMesh::DenseMatrix< T >::vector_mult_add().

390 {
391  libmesh_assert_equal_to (this->size(), vec.size());
392 
393  for (unsigned int i=0; i<this->size(); i++)
394  (*this)(i) += factor*vec(i);
395 }
template<typename T >
template<typename T2 >
CompareTypes< T, T2 >::supertype libMesh::DenseVector< T >::dot ( const DenseVector< T2 > &  vec) const
inline

Evaluate dot product with vec. In the complex-valued case, use the complex conjugate of vec.

Definition at line 400 of file dense_vector.h.

References libMesh::libmesh_conj().

401 {
402  libmesh_assert_equal_to (this->size(), vec.size());
403 
404  typename CompareTypes<T, T2>::supertype val = 0.;
405 
406  for (unsigned int i=0; i<this->size(); i++)
407  val += (*this)(i)*libmesh_conj(vec(i));
408 
409  return val;
410 }
template<typename T>
virtual T libMesh::DenseVector< T >::el ( const unsigned int  i) const
inlinevirtual
Returns
the (i) element of the vector.

Implements libMesh::DenseVectorBase< T >.

Definition at line 103 of file dense_vector.h.

103 { return (*this)(i); }
template<typename T>
virtual T& libMesh::DenseVector< T >::el ( const unsigned int  i)
inlinevirtual
Returns
the (i) element of the vector as a writeable reference.

Implements libMesh::DenseVectorBase< T >.

Definition at line 108 of file dense_vector.h.

108 { return (*this)(i); }
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 569 of file dense_vector.h.

571 {
572  libmesh_assert_less_equal ( sub_n, this->size() );
573 
574  dest.resize(sub_n);
575  for(unsigned int i=0; i<sub_n; i++)
576  dest(i) = (*this)(i);
577 }
template<typename T>
std::vector<T>& libMesh::DenseVector< T >::get_values ( )
inline

Access to the values array. This should be used with caution but can be used to speed up code compilation significantly.

Definition at line 230 of file dense_vector.h.

Referenced by libMesh::DenseMatrix< T >::_lu_back_substitute_lapack(), libMesh::DenseMatrix< T >::_matvec_blas(), libMesh::EpetraVector< T >::add_vector(), and libMesh::FEMContext::pre_fe_reinit().

230 { return _val; }
template<typename T>
const std::vector<T>& libMesh::DenseVector< T >::get_values ( ) const
inline

Access to the values array. This should be used with caution but can be used to speed up code compilation significantly.

Definition at line 237 of file dense_vector.h.

237 { return _val; }
template<typename T >
template<typename T2 >
CompareTypes< T, T2 >::supertype libMesh::DenseVector< T >::indefinite_dot ( const DenseVector< T2 > &  vec) const
inline

Evaluate dot product with vec. In the complex-valued case, do not use the complex conjugate of vec.

Definition at line 415 of file dense_vector.h.

416 {
417  libmesh_assert_equal_to (this->size(), vec.size());
418 
419  typename CompareTypes<T, T2>::supertype val = 0.;
420 
421  for (unsigned int i=0; i<this->size(); i++)
422  val += (*this)(i)*(vec(i));
423 
424  return val;
425 }
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.

Definition at line 525 of file dense_vector.h.

References std::abs(), and libMesh::Real.

526 {
527  Real my_norm = 0.;
528  for (unsigned int i=0; i!=this->size(); i++)
529  {
530  my_norm += std::abs((*this)(i));
531  }
532  return my_norm;
533 }
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 elements.

Definition at line 539 of file dense_vector.h.

References libMesh::TensorTools::norm_sq(), and libMesh::Real.

540 {
541  Real my_norm = 0.;
542  for (unsigned int i=0; i!=this->size(); i++)
543  {
544  my_norm += TensorTools::norm_sq((*this)(i));
545  }
546  return sqrt(my_norm);
547 }
template<typename T >
Real libMesh::DenseVector< T >::linfty_norm ( ) const
inline
Returns
the maximum absolute value of the elements of this vector, which is the $l_\infty$-norm of a vector.

Definition at line 553 of file dense_vector.h.

References libMesh::TensorTools::norm_sq(), and libMesh::Real.

554 {
555  if (!this->size())
556  return 0.;
557  Real my_norm = TensorTools::norm_sq((*this)(0));
558 
559  for (unsigned int i=1; i!=this->size(); i++)
560  {
561  Real current = TensorTools::norm_sq((*this)(i));
562  my_norm = (my_norm > current? my_norm : current);
563  }
564  return sqrt(my_norm);
565 }
template<typename T >
Real libMesh::DenseVector< T >::max ( ) const
inline
Returns
the maximum element in the vector. In case of complex numbers, this returns the maximum Real part.

Definition at line 508 of file dense_vector.h.

References libMesh::libmesh_assert(), libMesh::libmesh_real(), and libMesh::Real.

509 {
510  libmesh_assert (this->size());
511  Real my_max = libmesh_real((*this)(0));
512 
513  for (unsigned int i=1; i!=this->size(); i++)
514  {
515  Real current = libmesh_real((*this)(i));
516  my_max = (my_max > current? my_max : current);
517  }
518  return my_max;
519 }
template<typename T >
Real libMesh::DenseVector< T >::min ( ) const
inline
Returns
the minimum element in the vector. In case of complex numbers, this returns the minimum Real part.

Definition at line 491 of file dense_vector.h.

References libMesh::libmesh_assert(), libMesh::libmesh_real(), and libMesh::Real.

492 {
493  libmesh_assert (this->size());
494  Real my_min = libmesh_real((*this)(0));
495 
496  for (unsigned int i=1; i!=this->size(); i++)
497  {
498  Real current = libmesh_real((*this)(i));
499  my_min = (my_min < current? my_min : current);
500  }
501  return my_min;
502 }
template<typename T >
template<typename T2 >
bool libMesh::DenseVector< T >::operator!= ( const DenseVector< T2 > &  vec) const
inline

Tests if vec is not exactly equal to this vector.

Definition at line 446 of file dense_vector.h.

447 {
448  libmesh_assert_equal_to (this->size(), vec.size());
449 
450  for (unsigned int i=0; i<this->size(); i++)
451  if ((*this)(i) != vec(i))
452  return true;
453 
454  return false;
455 }
template<typename T >
T libMesh::DenseVector< T >::operator() ( const unsigned int  i) const
inline
Returns
the (i) element of the vector.

Definition at line 343 of file dense_vector.h.

344 {
345  libmesh_assert_less (i, _val.size());
346 
347  return _val[i];
348 }
template<typename T >
T & libMesh::DenseVector< T >::operator() ( const unsigned int  i)
inline
Returns
the (i,j) element of the vector as a writeable reference.

Definition at line 354 of file dense_vector.h.

355 {
356  libmesh_assert_less (i, _val.size());
357 
358  return _val[i];
359 }
template<typename T>
DenseVector< T > & libMesh::DenseVector< T >::operator*= ( const T  factor)
inline

Multiplies every element in the vector by factor.

Definition at line 375 of file dense_vector.h.

References libMesh::MeshTools::Modification::scale().

376 {
377  this->scale(factor);
378  return *this;
379 }
template<typename T >
template<typename T2 >
DenseVector< T > & libMesh::DenseVector< T >::operator+= ( const DenseVector< T2 > &  vec)
inline

Adds vec to this vector.

Definition at line 462 of file dense_vector.h.

463 {
464  libmesh_assert_equal_to (this->size(), vec.size());
465 
466  for (unsigned int i=0; i<this->size(); i++)
467  (*this)(i) += vec(i);
468 
469  return *this;
470 }
template<typename T >
template<typename T2 >
DenseVector< T > & libMesh::DenseVector< T >::operator-= ( const DenseVector< T2 > &  vec)
inline

Subtracts vec from this vector.

Definition at line 477 of file dense_vector.h.

478 {
479  libmesh_assert_equal_to (this->size(), vec.size());
480 
481  for (unsigned int i=0; i<this->size(); i++)
482  (*this)(i) -= vec(i);
483 
484  return *this;
485 }
template<typename T >
template<typename T2 >
DenseVector< T > & libMesh::DenseVector< T >::operator= ( const DenseVector< T2 > &  other_vector)
inline

Assignment operator.

Definition at line 293 of file dense_vector.h.

294 {
295  // _val = other_vector._val;
296 
297  const std::vector<T2> &other_vals = other_vector.get_values();
298 
299  _val.clear();
300  _val.reserve(other_vals.size());
301 
302  for (unsigned int i=0; i<other_vals.size(); i++)
303  _val.push_back(other_vals[i]);
304 
305  return *this;
306 }
template<typename T >
template<typename T2 >
bool libMesh::DenseVector< T >::operator== ( const DenseVector< T2 > &  vec) const
inline

Tests if vec is exactly equal to this vector.

Definition at line 430 of file dense_vector.h.

431 {
432  libmesh_assert_equal_to (this->size(), vec.size());
433 
434  for (unsigned int i=0; i<this->size(); i++)
435  if ((*this)(i) != vec(i))
436  return false;
437 
438  return true;
439 }
template<typename T >
void libMesh::DenseVectorBase< T >::print ( std::ostream &  os) const
inherited

Pretty-print the vector to stdout.

Definition at line 62 of file dense_vector_base.C.

63 {
64  for (unsigned int i=0; i<this->size(); i++)
65  os << std::setw(8)
66  << this->el(i)
67  << std::endl;
68 }
template<typename T >
void libMesh::DenseVectorBase< T >::print_scientific ( std::ostream &  os) const
inherited

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

Definition at line 30 of file dense_vector_base.C.

31 {
32 #ifndef LIBMESH_BROKEN_IOSTREAM
33 
34  // save the initial format flags
35  std::ios_base::fmtflags os_flags = os.flags();
36 
37  // Print the vector entries.
38  for (unsigned int i=0; i<this->size(); i++)
39  os << std::setw(10)
40  << std::scientific
41  << std::setprecision(8)
42  << this->el(i)
43  << std::endl;
44 
45  // reset the original format flags
46  os.flags(os_flags);
47 
48 #else
49 
50  // Print the matrix entries.
51  for (unsigned int i=0; i<this->size(); i++)
52  os << std::setprecision(8)
53  << this->el(i)
54  << std::endl;
55 
56 #endif
57 }
template<typename T>
void libMesh::DenseVector< T >::scale ( const T  factor)
inline

Multiplies every element in the vector by factor.

Definition at line 365 of file dense_vector.h.

366 {
367  for (unsigned int i=0; i<_val.size(); i++)
368  _val[i] *= factor;
369 }
template<typename T>
virtual unsigned int libMesh::DenseVector< T >::size ( ) const
inlinevirtual
Returns
the size of the vector.

Implements libMesh::DenseVectorBase< T >.

Definition at line 81 of file dense_vector.h.

Referenced by libMesh::DenseMatrix< T >::_lu_back_substitute(), libMesh::DenseMatrix< T >::_matvec_blas(), libMesh::HPCoarsenTest::add_projection(), libMesh::DistributedVector< T >::add_vector(), libMesh::LaspackVector< T >::add_vector(), libMesh::EigenSparseVector< T >::add_vector(), libMesh::EpetraVector< T >::add_vector(), libMesh::PetscVector< T >::add_vector(), libMesh::DiscontinuityMeasure::boundary_side_integration(), libMesh::KellyErrorEstimator::boundary_side_integration(), libMesh::DofMap::constrain_element_dyad_matrix(), libMesh::DofMap::constrain_element_matrix_and_vector(), libMesh::DofMap::constrain_element_vector(), libMesh::EulerSolver::element_residual(), libMesh::Euler2Solver::element_residual(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::DofMap::heterogenously_constrain_element_matrix_and_vector(), libMesh::DofMap::heterogenously_constrain_element_vector(), libMesh::LaspackVector< T >::insert(), libMesh::DistributedVector< T >::insert(), libMesh::EigenSparseVector< T >::insert(), libMesh::EpetraVector< T >::insert(), libMesh::PetscVector< T >::insert(), libMesh::LaplacianErrorEstimator::internal_side_integration(), libMesh::DiscontinuityMeasure::internal_side_integration(), libMesh::KellyErrorEstimator::internal_side_integration(), libMesh::ConstFunction< Output >::operator()(), libMesh::ConstFEMFunction< Output >::operator()(), libMesh::WrappedFunction< Output >::operator()(), libMesh::ParsedFunction< Output >::operator()(), libMesh::HPCoarsenTest::select_refinement(), libMesh::EulerSolver::side_residual(), libMesh::Euler2Solver::side_residual(), libMesh::DenseMatrix< T >::vector_mult(), libMesh::DenseMatrix< T >::vector_mult_add(), and libMesh::DenseMatrix< T >::vector_mult_transpose().

81  {
82  return libmesh_cast_int<unsigned int>(_val.size());
83  }
template<typename T>
void libMesh::DenseVector< T >::swap ( DenseVector< T > &  other_vector)
inline

STL-like swap method

Definition at line 312 of file dense_vector.h.

Referenced by libMesh::EulerSolver::element_residual(), libMesh::Euler2Solver::element_residual(), libMesh::EulerSolver::side_residual(), and libMesh::Euler2Solver::side_residual().

313 {
314  _val.swap(other_vector._val);
315 }

Member Data Documentation

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

The actual data values, stored as a 1D array.

Definition at line 244 of file dense_vector.h.

Referenced by libMesh::DenseVector< T >::DenseVector(), libMesh::DenseVector< Number >::get_values(), and libMesh::DenseVector< Number >::size().


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

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:58 UTC

Hosted By:
SourceForge.net Logo