Tpetra parallel linear algebra  Version of the Day
Public Member Functions | List of all members
Tpetra::Experimental::LittleVector< Scalar, LO > Class Template Reference

Nonowning view of a set of degrees of freedom corresponding to a mesh point in a block vector or multivector. More...

#include <Tpetra_Experimental_BlockView.hpp>

Public Member Functions

 LittleVector (Scalar *const A, const LO blockSize, const LO stride)
 Constructor. More...
 
template<class T >
 LittleVector (T *const A, const LO blockSize, const LO stride, typename std::enable_if< !std::is_same< Scalar, T >::value &&std::is_convertible< Scalar, T >::value &&sizeof(Scalar)==sizeof(T), int * >::type ignoreMe=NULL)
 Constructor that takes an impl_scalar_type pointer. More...
 
Scalar * getRawPtr () const
 Pointer to the block's entries. More...
 
LO getBlockSize () const
 The block size (number of degrees of freedom per mesh point). More...
 
LO getStride () const
 Stride between consecutive entries. More...
 
impl_scalar_type & operator() (const LO i) const
 Reference to entry (i) of the vector. More...
 
template<class LittleVectorType >
void update (const Scalar &alpha, const LittleVectorType &X) const
 *this := *this + alpha * X. More...
 
template<class LittleVectorType >
void assign (const LittleVectorType &X) const
 *this := X. More...
 
void scale (const Scalar &alpha) const
 (*this)(i,j) := alpha * (*this)(i,j) for all (i,j). More...
 
void fill (const Scalar &alpha) const
 (*this)(i,j) := alpha for all (i,j). More...
 
template<class LittleVectorType >
void absmax (const LittleVectorType &X) const
 (*this)(i,j) := max(abs((*this)(i,j)), abs(X(i,j))) for all (i,j). More...
 
template<class LittleVectorType >
bool equal (const LittleVectorType &X) const
 true if and only if all entries of this equal all entries of X. More...
 
template<class LittleBlockType , class LittleVectorType >
void matvecUpdate (const Scalar &alpha, const LittleBlockType &A, const LittleVectorType &X) const
 (*this) := (*this) + alpha * A * X (matrix-vector multiply). More...
 

Detailed Description

template<class Scalar, class LO>
class Tpetra::Experimental::LittleVector< Scalar, LO >

Nonowning view of a set of degrees of freedom corresponding to a mesh point in a block vector or multivector.

Template Parameters
ScalarThe type of entries.
LOThe type of local indices. See the documentation of the first template parameter of Map for requirements.

"Little" means local (not distributed over multiple MPI processes; stored to maximize locality) and small (think length 3, not length 1000).

The Scalar template parameter may be const or nonconst. This is one reason why instance methods below that take a LittleVector accept it as a template parameter: that lets you add a const LittleVector (e.g., LittleVector<const double, int>) to a nonconst LittleVector (e.g., LittleVector<double, int>).

Definition at line 325 of file Tpetra_Experimental_BlockView.hpp.

Constructor & Destructor Documentation

template<class Scalar , class LO >
Tpetra::Experimental::LittleVector< Scalar, LO >::LittleVector ( Scalar *const  A,
const LO  blockSize,
const LO  stride 
)
inline

Constructor.

Parameters
A[in] Pointer to the vector's entries
blockSize[in] Dimension of the vector
stride[in] Stride between consecutive entries

Definition at line 338 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
template<class T >
Tpetra::Experimental::LittleVector< Scalar, LO >::LittleVector ( T *const  A,
const LO  blockSize,
const LO  stride,
typename std::enable_if< !std::is_same< Scalar, T >::value &&std::is_convertible< Scalar, T >::value &&sizeof(Scalar)==sizeof(T), int * >::type  ignoreMe = NULL 
)
inline

Constructor that takes an impl_scalar_type pointer.

Parameters
A[in] Pointer to the vector's entries, as impl_scalar_type* rather than Scalar*
blockSize[in] Dimension of the vector
stride[in] Stride between consecutive entries

While this constructor is templated on a type T, the intent is that T == impl_scalar_type. (We must template on T rather than using impl_scalar_type directly, because of how std::enable_if works.) The long, complicated std::enable_if expression ensures that this constructor only exists if Scalar differs from impl_scalar_type, but the two types are mutually compatible and have the same size. (They must be bitwise compatible, so that reinterpret_cast makes sense between them.)

Definition at line 361 of file Tpetra_Experimental_BlockView.hpp.

Member Function Documentation

template<class Scalar , class LO >
Scalar* Tpetra::Experimental::LittleVector< Scalar, LO >::getRawPtr ( ) const
inline

Pointer to the block's entries.

Definition at line 375 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
LO Tpetra::Experimental::LittleVector< Scalar, LO >::getBlockSize ( ) const
inline

The block size (number of degrees of freedom per mesh point).

Definition at line 380 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
LO Tpetra::Experimental::LittleVector< Scalar, LO >::getStride ( ) const
inline

Stride between consecutive entries.

Definition at line 385 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
impl_scalar_type& Tpetra::Experimental::LittleVector< Scalar, LO >::operator() ( const LO  i) const
inline

Reference to entry (i) of the vector.

Note
To Tpetra developers: This is returned as impl_scalar_type and not as Scalar, in order to avoid a lot of reinterpret_cast calls in the inner loop of the sparse matrix-vector multiply kernel of Tpetra::Experimental::BlockCrsMatrix. Any pair of types impl_scalar_type, Scalar used here should always be convertible in either direction, so the return type should not pose any issues in practice.

Definition at line 399 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
template<class LittleVectorType >
void Tpetra::Experimental::LittleVector< Scalar, LO >::update ( const Scalar &  alpha,
const LittleVectorType &  X 
) const
inline

*this := *this + alpha * X.

Definition at line 405 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
template<class LittleVectorType >
void Tpetra::Experimental::LittleVector< Scalar, LO >::assign ( const LittleVectorType &  X) const
inline

*this := X.

Definition at line 414 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
void Tpetra::Experimental::LittleVector< Scalar, LO >::scale ( const Scalar &  alpha) const
inline

(*this)(i,j) := alpha * (*this)(i,j) for all (i,j).

Definition at line 421 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
void Tpetra::Experimental::LittleVector< Scalar, LO >::fill ( const Scalar &  alpha) const
inline

(*this)(i,j) := alpha for all (i,j).

Definition at line 429 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
template<class LittleVectorType >
void Tpetra::Experimental::LittleVector< Scalar, LO >::absmax ( const LittleVectorType &  X) const
inline

(*this)(i,j) := max(abs((*this)(i,j)), abs(X(i,j))) for all (i,j).

Tpetra uses this operation to implement the ABSMAX CombineMode.

Definition at line 441 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
template<class LittleVectorType >
bool Tpetra::Experimental::LittleVector< Scalar, LO >::equal ( const LittleVectorType &  X) const
inline

true if and only if all entries of this equal all entries of X.

Definition at line 450 of file Tpetra_Experimental_BlockView.hpp.

template<class Scalar , class LO >
template<class LittleBlockType , class LittleVectorType >
void Tpetra::Experimental::LittleVector< Scalar, LO >::matvecUpdate ( const Scalar &  alpha,
const LittleBlockType &  A,
const LittleVectorType &  X 
) const
inline

(*this) := (*this) + alpha * A * X (matrix-vector multiply).

Definition at line 465 of file Tpetra_Experimental_BlockView.hpp.


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