42 #ifndef THYRA_DEFAULT_PRODUCT_VECTOR_SPACE_HPP 43 #define THYRA_DEFAULT_PRODUCT_VECTOR_SPACE_HPP 46 #include "Thyra_DefaultProductVectorSpace_decl.hpp" 47 #include "Thyra_DefaultProductVector.hpp" 48 #include "Thyra_DefaultProductMultiVector.hpp" 49 #include "Thyra_ProductMultiVectorBase.hpp" 50 #include "Teuchos_Workspace.hpp" 51 #include "Teuchos_dyn_cast.hpp" 60 template<
class Scalar>
62 : numBlocks_(-1), dim_(-1), isInCore_(false)
66 template<
class Scalar>
70 : numBlocks_(-1), dim_(-1)
76 template<
class Scalar>
85 const int nBlocks = vecSpaces_in.size();
87 TEUCHOS_TEST_FOR_EXCEPT( nBlocks == 0 );
89 bool overallHasInCoreView =
true;
90 for (
int k = 0; k < nBlocks; ++k) {
92 TEUCHOS_TEST_FOR_EXCEPTION(
93 vecSpaces_in[k].
get() == NULL, std::invalid_argument
94 ,
"Error, the smart pointer vecSpaces["<<k<<
"] can not be NULL!" 97 if (!vecSpaces_in[k]->hasInCoreView()) overallHasInCoreView =
false;
103 numBlocks_ = nBlocks;
104 vecSpaces_ = Teuchos::rcp(
new vecSpaces_t);
105 *vecSpaces_ = vecSpaces_in;
106 vecSpacesOffsets_ = Teuchos::rcp(
new vecSpacesOffsets_t(nBlocks+1));
107 (*vecSpacesOffsets_)[0] = 0;
109 for(
int k = 1; k <= nBlocks; ++k ) {
110 const Ordinal dim_km1 = vecSpaces_in[k-1]->dim();
111 (*vecSpacesOffsets_)[k] = (*vecSpacesOffsets_)[k-1] + dim_km1;
114 isInCore_ = overallHasInCoreView;
119 template<
class Scalar>
124 TEUCHOS_TEST_FOR_EXCEPT(!is_null(vecSpaces_in));
125 vecSpaces_ = Teuchos::null;
126 vecSpacesOffsets_ = Teuchos::null;
133 template<
class Scalar>
140 TEUCHOS_TEST_FOR_EXCEPTION(
141 !(0 <= i && i < this->dim()), std::out_of_range
142 ,
"VectorSpaceBlocked::get_vector_space_position(...): Error, i = " 143 << i <<
" is not in range [0,"<<(this->dim()-1)<<
"]" 146 *kth_vector_space = 0;
147 *kth_global_offset = 0;
148 while( *kth_vector_space < numBlocks_ ) {
149 const Ordinal off_kp1 = (*vecSpacesOffsets_)[*kth_vector_space+1];
151 *kth_global_offset = (*vecSpacesOffsets_)[*kth_vector_space];
154 ++(*kth_vector_space);
156 TEUCHOS_TEST_FOR_EXCEPT( !(*kth_vector_space < numBlocks_) );
163 template<
class Scalar>
170 template<
class Scalar>
171 Teuchos::RCP<const VectorSpaceBase<Scalar> >
174 TEUCHOS_TEST_FOR_EXCEPT( k < 0 || numBlocks_ < k );
175 return (*vecSpaces_)[k];
182 template<
class Scalar>
189 template<
class Scalar>
194 using Teuchos::ptrFromRef;
195 using Teuchos::ptr_dynamic_cast;
197 const int nBlocks = this->numBlocks();
200 const Ptr<const ProductVectorSpaceBase<Scalar> > pvsb =
205 if( nBlocks != pvsb->numBlocks() )
207 for(
int i = 0; i < nBlocks; ++i ) {
208 if( !this->getBlock(i)->isCompatible(*pvsb->getBlock(i)) )
226 template<
class Scalar>
227 Teuchos::RCP< VectorBase<Scalar> >
230 return defaultProductVector<Scalar>(Teuchos::rcpFromRef(*
this));
234 template<
class Scalar>
240 const int nBlocks = this->numBlocks();
245 TEUCHOS_TEST_FOR_EXCEPT(
247 || nBlocks!=y.productSpace()->numBlocks()
250 Scalar scalarProd_rtn = Teuchos::ScalarTraits<Scalar>::zero();
251 for(
int k = 0; k < nBlocks; ++k )
252 scalarProd_rtn += (*vecSpaces_)[k]->scalarProd(
255 return scalarProd_rtn;
259 template<
class Scalar>
263 const ArrayView<Scalar> &scalarProds_out
267 using Teuchos::Workspace;
271 TEUCHOS_TEST_FOR_EXCEPT(is_null(scalarProds_out));
273 TEUCHOS_ASSERT_EQUALITY( as<Ordinal>(scalarProds_out.size()),
277 scalarProds_out[0] = this->scalarProd(*X_in.
col(0),*Y_in.
col(0));
281 Teuchos::WorkspaceStore* wss = Teuchos::get_default_workspace_store().get();
282 const int nBlocks = this->numBlocks();
287 TEUCHOS_TEST_FOR_EXCEPT( nBlocks!=X.
productSpace()->numBlocks() || nBlocks!=Y.productSpace()->numBlocks() );
289 Workspace<Scalar> _scalarProds_out(wss, m,
false);
290 std::fill( scalarProds_out.begin(), scalarProds_out.end(),
291 ScalarTraits<Scalar>::zero() );
292 for(
int k = 0; k < nBlocks; ++k ) {
293 (*vecSpaces_)[k]->scalarProds(
295 for(
int j = 0; j < m; ++j )
296 scalarProds_out[j] += _scalarProds_out[j];
301 template<
class Scalar>
304 const Range1D rng = full_range(rng_in,0,dim_-1);
306 int kth_vector_space = -1;
308 this->getVecSpcPoss(rng.lbound(),&kth_vector_space,&kth_global_offset);
310 TEUCHOS_TEST_FOR_EXCEPT( !( 0 <= kth_vector_space && kth_vector_space <= numBlocks_ ) );
312 if( rng.lbound() + rng.size() <= kth_global_offset + (*vecSpaces_)[kth_vector_space]->dim() ) {
313 return (*vecSpaces_)[kth_vector_space]->hasInCoreView(rng_in-kth_global_offset,viewType,strideType);
337 template<
class Scalar>
338 Teuchos::RCP< const VectorSpaceFactoryBase<Scalar> >
342 return (*vecSpaces_)[0]->smallVecSpcFcty();
343 return Teuchos::null;
347 template<
class Scalar>
348 Teuchos::RCP< MultiVectorBase<Scalar> >
351 return defaultProductMultiVector<Scalar>(Teuchos::rcpFromRef(*
this),
356 template<
class Scalar>
357 Teuchos::RCP< const VectorSpaceBase<Scalar> >
363 Teuchos::RCP<DefaultProductVectorSpace<Scalar> >
364 pvs = productVectorSpace<Scalar>();
365 pvs->numBlocks_ = numBlocks_;
366 pvs->vecSpaces_ = vecSpaces_;
367 pvs->vecSpacesOffsets_ = vecSpacesOffsets_;
369 pvs->isInCore_ = isInCore_;
377 template<
class Scalar>
380 std::ostringstream oss;
382 << Teuchos::Describable::description() <<
"{" 384 <<
",numBlocks="<<numBlocks_
390 template<
class Scalar>
392 Teuchos::FancyOStream &out_arg
393 ,
const Teuchos::EVerbosityLevel verbLevel
396 using Teuchos::includesVerbLevel;
397 using Teuchos::OSTab;
398 RCP<FancyOStream> out = rcpFromRef(out_arg);
400 if (includesVerbLevel(verbLevel, Teuchos::VERB_LOW,
true)) {
401 *out << this->description() << std::endl;
403 if (includesVerbLevel(verbLevel, Teuchos::VERB_MEDIUM) && numBlocks_ > 0) {
405 *out <<
"Constituent vector spaces V[0], V[1], ... V[numBlocks-1]:\n";
407 for(
int k = 0; k < numBlocks_; ++k ) {
408 *out <<
"V["<<k<<
"] = " << Teuchos::describe(*(*vecSpaces_)[k],verbLevel);
417 #endif // THYRA_DEFAULT_PRODUCT_VECTOR_SPACE_HPP Base interface for product multi-vectors.
Base interface for product vectors.
bool hasInCoreView(const Range1D &rng, const EViewType viewType, const EStrideType strideType) const
Returns true if all of the constituent vector spaces return true.
virtual RCP< const VectorBase< Scalar > > getVectorBlock(const int k) const =0
Returns a non-persisting const view of the (zero-based) kth block vector.
The view is a direct view of data and no copies are made.
virtual Teuchos::RCP< const ProductVectorSpaceBase< Scalar > > productSpace() const =0
Returns the associated product vector space that represents the range.
RCP< VectorBase< Scalar > > createMember() const
Returns a DefaultProductVector object.
RCP< const VectorSpaceFactoryBase< Scalar > > smallVecSpcFcty() const
Returns getBlock(0)->smallVecSpcFcty().
Scalar scalarProd(const VectorBase< Scalar > &x, const VectorBase< Scalar > &y) const
Returns the sum of the scalar products of the constituent vectors.
Abstract interface for objects that represent a space for vectors.
EViewType
Determines if a view is a direct view of data or a detached copy of data.
virtual Teuchos::RCP< const MultiVectorBase< Scalar > > getMultiVectorBlock(const int k) const =0
Returns a non-persisting const view of the (zero-based) kth block multi-vector.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
Prints the details about the constituent vector spaces.
Teuchos::Ordinal Ordinal
Type for the dimension of a vector space. `*.
Interface for a collection of column vectors called a multi-vector.
void getVecSpcPoss(Ordinal i, int *kth_vector_space, Ordinal *kth_global_offset) const
Get the position of the vector space object and its offset into a composite vector that owns the ith ...
RCP< const VectorBase< Scalar > > col(Ordinal j) const
Calls colImpl().
virtual void uninitialize(const ArrayView< RCP< const VectorSpaceBase< Scalar > > > &vecSpaces=Teuchos::null)
Uninitialize.
virtual void initialize(const ArrayView< const RCP< const VectorSpaceBase< Scalar > > > &vecSpaces)
Initialize with a list of constituent vector spaces.
RCP< const VectorSpaceBase< Scalar > > getBlock(const int k) const
Abstract interface for finite-dimensional dense vectors.
Ordinal dim() const
Returns the summation of the constituent vector spaces.
bool isCompatible(const VectorSpaceBase< Scalar > &vecSpc) const
Returns true only if also a product vector space and all constituent vectors are compatible.
virtual RCP< const VectorSpaceBase< Scalar > > domain() const =0
Return a smart pointer for the domain space for this operator.
RCP< const VectorSpaceBase< Scalar > > clone() const
Clones the object as promised.
void scalarProdsImpl(const MultiVectorBase< Scalar > &X, const MultiVectorBase< Scalar > &Y, const ArrayView< Scalar > &scalarProds) const
Returns the sum of the scalar products of each of the columns of the constituent multi-vectors.
EStrideType
Determine if data is unit stride or non-unit stride.
std::string description() const
Prints just the name DefaultProductVectorSpace along with the overall dimension and the number of blo...
virtual bool isCompatible(const VectorSpaceBase< Scalar > &vecSpc) const =0
Compare the compatibility of two vector spaces.
virtual Ordinal dim() const =0
Return the dimension of the vector space.
RCP< MultiVectorBase< Scalar > > createMembers(int numMembers) const
Returns a DefaultProductMultiVector object.
DefaultProductVectorSpace()
Default construct to uninitialized.