42 #ifndef THYRA_MULTI_VECTOR_DEFAULT_BASE_DEF_HPP 43 #define THYRA_MULTI_VECTOR_DEFAULT_BASE_DEF_HPP 46 #include "Thyra_MultiVectorDefaultBase_decl.hpp" 47 #include "Thyra_LinearOpDefaultBase.hpp" 48 #include "Thyra_MultiVectorStdOps.hpp" 49 #include "Thyra_VectorSpaceFactoryBase.hpp" 50 #include "Thyra_VectorSpaceBase.hpp" 51 #include "Thyra_VectorBase.hpp" 52 #include "Thyra_AssertOp.hpp" 53 #include "Thyra_DefaultColumnwiseMultiVector.hpp" 54 #include "RTOpPack_TOpAssignScalar.hpp" 55 #include "Teuchos_Workspace.hpp" 56 #include "Teuchos_Assert.hpp" 57 #include "Teuchos_as.hpp" 66 template<
class Scalar>
67 RCP<MultiVectorBase<Scalar> >
71 &l_domain = *this->domain(),
72 &l_range = *this->range();
73 RCP<MultiVectorBase<Scalar> >
74 copy = createMembers(l_range,l_domain.
dim());
75 ::Thyra::assign<Scalar>(copy.ptr(), *
this);
85 template<
class Scalar>
88 using Teuchos::tuple;
using Teuchos::null;
90 Thyra::applyOp<Scalar>(assign_scalar_op,
91 ArrayView<Ptr<const MultiVectorBase<Scalar> > >(null),
96 template<
class Scalar>
97 RCP<const MultiVectorBase<Scalar> >
100 using Teuchos::Workspace;
102 Teuchos::WorkspaceStore *wss = Teuchos::get_default_workspace_store().get();
106 const Range1D colRng = Teuchos::full_range(colRng_in,0,dimDomain-1);
107 if( colRng.lbound() == 0 && as<Ordinal>(colRng.ubound()) == dimDomain-1 )
108 return Teuchos::rcp(
this,
false);
109 if( colRng.size() ) {
111 Workspace< RCP< VectorBase<Scalar> > > col_vecs(wss,colRng.size());
112 for(
Ordinal j = colRng.lbound(); j <= colRng.ubound(); ++j )
113 col_vecs[j-colRng.lbound()] = Teuchos::rcp_const_cast<
VectorBase<Scalar> >(this->col(j));
116 this->range(),l_range.
smallVecSpcFcty()->createVecSpc(colRng.size()),col_vecs
120 return Teuchos::null;
124 template<
class Scalar>
125 RCP<MultiVectorBase<Scalar> >
128 using Teuchos::Workspace;
130 Teuchos::WorkspaceStore *wss = Teuchos::get_default_workspace_store().get();
134 const Range1D colRng = Teuchos::full_range(colRng_in,0,dimDomain-1);
135 if( colRng.lbound() == 0 && as<Ordinal>(colRng.ubound()) == dimDomain-1 )
136 return Teuchos::rcp(
this,
false);
137 if( colRng.size() ) {
139 Workspace< RCP< VectorBase<Scalar> > > col_vecs(wss,colRng.size());
140 for(
Ordinal j = colRng.lbound(); j <= colRng.ubound(); ++j )
141 col_vecs[j-colRng.lbound()] = this->col(j);
144 this->range(),l_range.
smallVecSpcFcty()->createVecSpc(colRng.size()),col_vecs
148 return Teuchos::null;
152 template<
class Scalar>
153 RCP<const MultiVectorBase<Scalar> >
155 const ArrayView<const int> &cols
158 using Teuchos::Workspace;
159 Teuchos::WorkspaceStore *wss = Teuchos::get_default_workspace_store().get();
161 const int numCols = cols.size();
165 const char msg_err[] =
"MultiVectorDefaultBase<Scalar>::subView(numCols,cols[]): Error!";
166 TEUCHOS_TEST_FOR_EXCEPTION( numCols < 1 || dimDomain < numCols, std::invalid_argument, msg_err );
169 Workspace< RCP< VectorBase<Scalar> > > col_vecs(wss,numCols);
170 for(
int k = 0; k < numCols; ++k ) {
171 const int col_k = cols[k];
173 TEUCHOS_TEST_FOR_EXCEPTION(
174 !( 0 <= col_k && col_k < dimDomain ), std::invalid_argument
175 ,msg_err <<
" col["<<k<<
"] = " << col_k <<
" is not in the range [0,"<<(dimDomain-1)<<
"]!" 182 this->range(), l_range.
smallVecSpcFcty()->createVecSpc(numCols), col_vecs
188 template<
class Scalar>
189 RCP<MultiVectorBase<Scalar> >
191 const ArrayView<const int> &cols
194 using Teuchos::Workspace;
195 Teuchos::WorkspaceStore *wss = Teuchos::get_default_workspace_store().get();
197 const int numCols = cols.size();
201 const char msg_err[] =
"MultiVectorDefaultBase<Scalar>::subView(numCols,cols[]): Error!";
202 TEUCHOS_TEST_FOR_EXCEPTION( numCols < 1 || dimDomain < numCols, std::invalid_argument, msg_err );
205 Workspace< RCP< VectorBase<Scalar> > > col_vecs(wss,numCols);
206 for(
int k = 0; k < numCols; ++k ) {
207 const int col_k = cols[k];
209 TEUCHOS_TEST_FOR_EXCEPTION(
210 !( 0 <= col_k && col_k < dimDomain ), std::invalid_argument
211 ,msg_err <<
" col["<<k<<
"] = " << col_k <<
" is not in the range [0,"<<(dimDomain-1)<<
"]!" 214 col_vecs[k] = this->col(col_k);
218 this->range(), l_range.
smallVecSpcFcty()->createVecSpc(numCols), col_vecs
224 template<
class Scalar>
229 const ArrayView<
const Ptr<RTOpPack::ReductTarget> > &reduct_objs,
230 const Ordinal prim_global_offset_in
234 using Teuchos::Workspace;
236 Teuchos::WorkspaceStore* wss = Teuchos::get_default_workspace_store().get();
238 const int num_multi_vecs = multi_vecs.size();
239 const int num_targ_multi_vecs = targ_multi_vecs.size();
254 Workspace<RCP<const VectorBase<Scalar> > > vecs_s(wss, num_multi_vecs);
255 Workspace<Ptr<const VectorBase<Scalar> > > vecs(wss, num_multi_vecs);
256 Workspace<RCP<VectorBase<Scalar> > > targ_vecs_s(wss, num_targ_multi_vecs);
257 Workspace<Ptr<VectorBase<Scalar> > > targ_vecs(wss, num_targ_multi_vecs);
259 for(
Ordinal j = 0; j < sec_dim; ++j) {
261 {
for(
Ordinal k = 0; k < as<Ordinal>(num_multi_vecs); ++k) {
262 vecs_s[k] = multi_vecs[k]->col(j);
263 vecs[k] = vecs_s[k].ptr();
265 {
for(
Ordinal k = 0; k < as<Ordinal>(num_targ_multi_vecs); ++k) {
266 targ_vecs_s[k] = targ_multi_vecs[k]->col(j);
267 targ_vecs[k] = targ_vecs_s[k].ptr();
273 targ_vecs().getConst(),
274 reduct_objs.size() ? reduct_objs[j] : Ptr<RTOpPack::ReductTarget>(),
275 prim_global_offset_in);
283 template<
class Scalar>
289 const Ptr<RTOpPack::ReductTarget> &reduct_obj,
290 const Ordinal prim_global_offset_in
294 using Teuchos::Workspace;
295 Teuchos::WorkspaceStore* wss = Teuchos::get_default_workspace_store().get();
306 const int reduct_objs_size = (!is_null(reduct_obj) ? sec_dim : 0);
307 Workspace<RCP<RTOpPack::ReductTarget> > rcp_reduct_objs(wss, reduct_objs_size);
308 Workspace<Ptr<RTOpPack::ReductTarget> > reduct_objs(wss, reduct_objs_size);
309 if (!is_null(reduct_obj)) {
310 for(
Ordinal k = 0; k < sec_dim; ++k) {
312 reduct_objs[k] = rcp_reduct_objs[k].ptr();
318 prim_op, multi_vecs, targ_multi_vecs, reduct_objs,
319 prim_global_offset_in);
323 if (!is_null(reduct_obj)) {
324 for (
Ordinal k = 0; k < sec_dim; ++k) {
331 template<
class Scalar>
339 rangeDim = this->range()->dim(),
340 domainDim = this->domain()->dim();
342 rowRng = rowRng_in.full_range() ?
Range1D(0,rangeDim-1) : rowRng_in,
343 colRng = colRng_in.full_range() ?
Range1D(0,domainDim-1) : colRng_in;
345 TEUCHOS_TEST_FOR_EXCEPTION(
346 !(rowRng.ubound() < rangeDim), std::out_of_range
347 ,
"MultiVectorDefaultBase<Scalar>::acquireDetachedMultiVectorViewImpl(...): Error, rowRng = [" 348 <<rowRng.lbound()<<
","<<rowRng.ubound()<<
"] is not in the range = [0," 351 TEUCHOS_TEST_FOR_EXCEPTION(
352 !(colRng.ubound() < domainDim), std::out_of_range
353 ,
"MultiVectorDefaultBase<Scalar>::acquireDetachedMultiVectorViewImpl(...): Error, colRng = [" 354 <<colRng.lbound()<<
","<<colRng.ubound()<<
"] is not in the range = [0," 355 <<(domainDim-1)<<
"]!" 359 const ArrayRCP<Scalar> values = Teuchos::arcp<Scalar>(rowRng.size() * colRng.size());
362 for(
int k = colRng.lbound(); k <= colRng.ubound(); ++k ) {
363 RCP<const VectorBase<Scalar> > col_k = this->col(k);
364 col_k->acquireDetachedView( rowRng, &sv );
365 for(
int i = 0; i < rowRng.size(); ++i )
366 values[ i + k*rowRng.size() ] = sv[i];
367 col_k->releaseDetachedView( &sv );
381 template<
class Scalar>
391 template<
class Scalar>
412 template<
class Scalar>
418 TEUCHOS_TEST_FOR_EXCEPTION(
419 sub_mv==NULL, std::logic_error,
420 "MultiVectorDefaultBase<Scalar>::commitNonconstDetachedMultiVectorViewImpl(...): Error!" 427 RCP<VectorBase<Scalar> > col_k = this->col(k);
428 col_k->acquireDetachedView( rowRng, &msv );
429 for(
int i = 0; i < rowRng.size(); ++i )
430 msv[i] = sub_mv->
values()[ i + k*rowRng.size() ];
431 col_k->commitDetachedView( &msv );
441 #endif // THYRA_MULTI_VECTOR_DEFAULT_BASE_DEF_HPP RCP< MultiVectorBase< Scalar > > nonconstContigSubViewImpl(const Range1D &colRng)
virtual void mvMultiReductApplyOpImpl(const RTOpPack::RTOpT< Scalar > &primary_op, const ArrayView< const Ptr< const MultiVectorBase< Scalar > > > &multi_vecs, const ArrayView< const Ptr< MultiVectorBase< Scalar > > > &targ_multi_vecs, const ArrayView< const Ptr< RTOpPack::ReductTarget > > &reduct_objs, const Ordinal primary_global_offset) const
RCP< const MultiVectorBase< Scalar > > contigSubViewImpl(const Range1D &colRng) const
virtual void acquireNonconstDetachedMultiVectorViewImpl(const Range1D &rowRng, const Range1D &colRng, RTOpPack::SubMultiVectorView< Scalar > *sub_mv)
virtual void mvSingleReductApplyOpImpl(const RTOpPack::RTOpT< Scalar > &primary_op, const RTOpPack::RTOpT< Scalar > &secondary_op, const ArrayView< const Ptr< const MultiVectorBase< Scalar > > > &multi_vecs, const ArrayView< const Ptr< MultiVectorBase< Scalar > > > &targ_multi_vecs, const Ptr< RTOpPack::ReductTarget > &reduct_obj, const Ordinal primary_global_offset) const
Abstract interface for objects that represent a space for vectors.
void assignImpl(Scalar alpha)
virtual void releaseDetachedMultiVectorViewImpl(RTOpPack::ConstSubMultiVectorView< Scalar > *sub_mv) const
Teuchos::Ordinal Ordinal
Type for the dimension of a vector space. `*.
Interface for a collection of column vectors called a multi-vector.
virtual void acquireDetachedMultiVectorViewImpl(const Range1D &rowRng, const Range1D &colRng, RTOpPack::ConstSubMultiVectorView< Scalar > *sub_mv) const
Teuchos::RCP< ReductTarget > reduct_obj_create() const
virtual RCP< MultiVectorBase< Scalar > > clone_mv() const
virtual RCP< const VectorSpaceFactoryBase< Scalar > > smallVecSpcFcty() const =0
Return a VectorSpaceFactoryBase object for the creation of (usually serial) vector spaces with a smal...
Abstract interface for finite-dimensional dense vectors.
Default subclass for MultiVectorBase implemented using columns of separate abstract vectors...
virtual void commitNonconstDetachedMultiVectorViewImpl(RTOpPack::SubMultiVectorView< Scalar > *sub_mv)
RCP< MultiVectorBase< Scalar > > nonconstNonContigSubViewImpl(const ArrayView< const int > &cols)
Ordinal numSubCols() const
const ArrayRCP< Scalar > values() const
void initialize(Ordinal globalOffset_in, Ordinal subDim_in, Ordinal colOffset_in, Ordinal numSubCols_in, const ArrayRCP< const Scalar > &values_in, Ordinal leadingDim_in)
Ordinal globalOffset() const
RCP< const MultiVectorBase< Scalar > > nonContigSubViewImpl(const ArrayView< const int > &cols) const
void reduce_reduct_objs(const ReductTarget &in_reduct_obj, const Ptr< ReductTarget > &inout_reduct_obj) const
Ordinal colOffset() const
virtual Ordinal dim() const =0
Return the dimension of the vector space.