AbstractLinAlgPack: C++ Interfaces For Vectors, Matrices And Related Linear Algebra Objects  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Public Types | Public Member Functions | List of all members
AbstractLinAlgPack::VectorMutableBlocked Class Reference

Concrete subclass for a blocked vector. More...

#include <AbstractLinAlgPack_VectorMutableBlocked.hpp>

Inheritance diagram for AbstractLinAlgPack::VectorMutableBlocked:
Inheritance graph
[legend]

Public Types

typedef Teuchos::RCP< const
VectorSpaceBlocked
vec_space_comp_ptr_t
 
- Public Types inherited from AbstractLinAlgPack::Vector
typedef Teuchos::RCP< const
Vector
vec_ptr_t
 
typedef Teuchos::RCP
< VectorMutable
vec_mut_ptr_t
 

Public Member Functions

 VectorMutableBlocked (VectorMutable::vec_mut_ptr_t *vecs, const vec_space_comp_ptr_t &vec_space)
 Calls this->initialize(). More...
 
void initialize (VectorMutable::vec_mut_ptr_t *vecs, const vec_space_comp_ptr_t &vec_space)
 Initialize given a list of constituent vectors. More...
 
const VectorSpaceBlockedblock_space () const
 Return the underlying block vector space. More...
 
const Vectorget_vector (int k) const
 Get the kth (zero based) constituent vector. More...
 
VectorMutableget_vector (int k)
 Get the kth (zero based) constituent vector. More...
 
- Public Member Functions inherited from AbstractLinAlgPack::VectorMutable
virtual VectorMutableoperator= (const VectorMutable &v)
 Default implementation calls operator=((const &Vector)v). More...
 
vec_mut_ptr_t sub_view (const index_type &l, const index_type &u)
 Inline member function that simply calls this->sub_view(Range1D(l,u)). More...
 
virtual void zero ()
 Zeros this vector. More...
 
virtual void get_sub_vector (const Range1D &rng, RTOpPack::MutableSubVector *sub_vec)
 Get a mutable explicit view of a sub-vector. More...
 
virtual void commit_sub_vector (RTOpPack::MutableSubVector *sub_vec)
 Free a mutable explicit view of a sub-vector. More...
 
virtual void Vp_StMtV (value_type alpha, const GenPermMatrixSlice &P, BLAS_Cpp::Transp P_trans, const Vector &x, value_type beta)
 Perform a gather or scatter operation with a vector. More...
 
vec_ptr_t sub_view (const Range1D &rng) const
 Default implementation calls this->sub_view() (non-const) and then performs an cast to vec_ptr_t. More...
 
- Public Member Functions inherited from AbstractLinAlgPack::Vector
 Vector ()
 
virtual ~Vector ()
 
vec_ptr_t sub_view (const index_type &l, const index_type &u) const
 Inline member function that simply calls this->sub_view(Range1D(l,u)). More...
 
virtual vec_mut_ptr_t clone () const
 Create a clone of this vector objet. More...
 
virtual value_type norm_2 () const
 Two norm. ||v||_2 = sqrt( sum( v(i)^2, i = 1,,,this->dim() ) ) More...
 

Overridden form Vector

index_type dim () const
 
const VectorSpacespace () const
 
void apply_op (const RTOpPack::RTOp &op, const size_t num_vecs, const Vector *vecs[], const size_t num_targ_vecs, VectorMutable *targ_vecs[], RTOpPack::ReductTarget *reduct_obj, const index_type first_ele, const index_type sub_dim, const index_type global_offset) const
 
index_type nz () const
 
std::ostream & output (std::ostream &out, bool print_dim, bool newline, index_type global_offset) const
 
value_type get_ele (index_type i) const
 
value_type norm_1 () const
 
value_type norm_inf () const
 
value_type inner_product (const Vector &) const
 
void get_sub_vector (const Range1D &rng, RTOpPack::SubVector *sub_vec) const
 
void free_sub_vector (RTOpPack::SubVector *sub_vec) const
 
void has_changed () const
 

Overridden from VectorMutable

vec_mut_ptr_t sub_view (const Range1D &rng)
 
void axpy (value_type alpha, const Vector &x)
 
VectorMutableoperator= (value_type)
 
VectorMutableoperator= (const Vector &)
 
void set_ele (index_type i, value_type val)
 
void set_sub_vector (const RTOpPack::SparseSubVector &sub_vec)
 

Detailed Description

Concrete subclass for a blocked vector.

This is a near optimal implementation for a block vector in many circumstances.

A blocked vector is simply the concatonation of two or more vectors to form a larger vector. Suppose that we are given a list of pointers to vector objects v[0],v[1],,,v[p]. Then this vector object would be represented as:

        [ *v[0] ]
        [ *v[1] ]
this =  [  .    ]
        [ *v[p] ]

This vector subclass trys to implement all of the vector methods as efficiently as possble. For example, if a client requests a sub-view that corresponds to a whole constituent sub-vector, then it will return the unadorned vector such as this->sub_view(1,v[0]->dim()).get() == v[0]. However, if some sub-view is requested that overlaps two or more constitient vectors, then there is no choice but to return a VectorMutableBlocked object from sub_view().

There are some situations where this vector subclass will be inefficient. For example, suppose that the constituent vectors v[0],v[1],,,v[p] are parallel vectors with elements owned by unique processes such that an reduction/transformation operator could be applied in parallel to all the constituent vectors in the same time (wall clock) that it takes to apply an operator to one constituent vector. In this case, this vector subclass could not take advantage of this parallelism and therefore, a more specialized "parallel" block vector class should be created (using threads for instance). Alternatively, the implementation of this vector subclass could be modified to take advantage of threads if they are available.

Definition at line 81 of file AbstractLinAlgPack_VectorMutableBlocked.hpp.

Member Typedef Documentation

Definition at line 86 of file AbstractLinAlgPack_VectorMutableBlocked.hpp.

Constructor & Destructor Documentation

AbstractLinAlgPack::VectorMutableBlocked::VectorMutableBlocked ( VectorMutable::vec_mut_ptr_t vecs,
const vec_space_comp_ptr_t vec_space 
)

Calls this->initialize().

Definition at line 66 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

Member Function Documentation

void AbstractLinAlgPack::VectorMutableBlocked::initialize ( VectorMutable::vec_mut_ptr_t vecs,
const vec_space_comp_ptr_t vec_space 
)

Initialize given a list of constituent vectors.

Preconditions:

  • vec_space.get() != NULL (throw std::logic_error).
  • vecs[k]->space().is_compatible(*vec_spaces->vector_spaces()[k]) == true, for k = 0...vec_space->num_vector_spaces().

Postconditions:

  • ToDo: Spell out for carefully!
Parameters
vecs[in] Array (size vec_space->num_vector_spaces()) of smart reference counted pointers to the constituent vectors.
vec_space[in] Block vector space object containing constituent vector spaces.

Definition at line 75 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

const VectorSpaceBlocked & AbstractLinAlgPack::VectorMutableBlocked::block_space ( ) const
inline

Return the underlying block vector space.

Definition at line 234 of file AbstractLinAlgPack_VectorMutableBlocked.hpp.

const Vector & AbstractLinAlgPack::VectorMutableBlocked::get_vector ( int  k) const
inline

Get the kth (zero based) constituent vector.

Definition at line 242 of file AbstractLinAlgPack_VectorMutableBlocked.hpp.

VectorMutable & AbstractLinAlgPack::VectorMutableBlocked::get_vector ( int  k)
inline

Get the kth (zero based) constituent vector.

Definition at line 250 of file AbstractLinAlgPack_VectorMutableBlocked.hpp.

index_type AbstractLinAlgPack::VectorMutableBlocked::dim ( ) const
virtual

Reimplemented from AbstractLinAlgPack::Vector.

Definition at line 93 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

const VectorSpace & AbstractLinAlgPack::VectorMutableBlocked::space ( ) const
virtual
void AbstractLinAlgPack::VectorMutableBlocked::apply_op ( const RTOpPack::RTOp &  op,
const size_t  num_vecs,
const Vector vecs[],
const size_t  num_targ_vecs,
VectorMutable targ_vecs[],
RTOpPack::ReductTarget reduct_obj,
const index_type  first_ele,
const index_type  sub_dim,
const index_type  global_offset 
) const
virtual
index_type AbstractLinAlgPack::VectorMutableBlocked::nz ( ) const
virtual

Reimplemented from AbstractLinAlgPack::Vector.

Definition at line 234 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

std::ostream & AbstractLinAlgPack::VectorMutableBlocked::output ( std::ostream &  out,
bool  print_dim,
bool  newline,
index_type  global_offset 
) const
virtual

Reimplemented from AbstractLinAlgPack::Vector.

Definition at line 245 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

value_type AbstractLinAlgPack::VectorMutableBlocked::get_ele ( index_type  i) const
virtual

Reimplemented from AbstractLinAlgPack::Vector.

Definition at line 265 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

value_type AbstractLinAlgPack::VectorMutableBlocked::norm_1 ( ) const
virtual

Reimplemented from AbstractLinAlgPack::Vector.

Definition at line 276 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

value_type AbstractLinAlgPack::VectorMutableBlocked::norm_inf ( ) const
virtual

Reimplemented from AbstractLinAlgPack::Vector.

Definition at line 291 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

value_type AbstractLinAlgPack::VectorMutableBlocked::inner_product ( const Vector v) const
virtual

Reimplemented from AbstractLinAlgPack::Vector.

Definition at line 306 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

void AbstractLinAlgPack::VectorMutableBlocked::get_sub_vector ( const Range1D rng,
RTOpPack::SubVector *  sub_vec 
) const
virtual

Reimplemented from AbstractLinAlgPack::Vector.

Definition at line 311 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

void AbstractLinAlgPack::VectorMutableBlocked::free_sub_vector ( RTOpPack::SubVector *  sub_vec) const
virtual

Reimplemented from AbstractLinAlgPack::Vector.

Definition at line 337 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.

VectorMutable::vec_mut_ptr_t AbstractLinAlgPack::VectorMutableBlocked::sub_view ( const Range1D rng)
virtual
void AbstractLinAlgPack::VectorMutableBlocked::axpy ( value_type  alpha,
const Vector x 
)
virtual
VectorMutable & AbstractLinAlgPack::VectorMutableBlocked::operator= ( value_type  alpha)
virtual
VectorMutable & AbstractLinAlgPack::VectorMutableBlocked::operator= ( const Vector v)
virtual
void AbstractLinAlgPack::VectorMutableBlocked::set_ele ( index_type  i,
value_type  val 
)
virtual
void AbstractLinAlgPack::VectorMutableBlocked::set_sub_vector ( const RTOpPack::SparseSubVector &  sub_vec)
virtual
void AbstractLinAlgPack::VectorMutableBlocked::has_changed ( ) const
protectedvirtual

Reimplemented from AbstractLinAlgPack::Vector.

Definition at line 358 of file AbstractLinAlgPack_VectorMutableBlocked.cpp.


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