52 #include "Teko_BlockedReordering.hpp" 55 #include "Teuchos_VerboseObject.hpp" 56 #include "Teuchos_RCP.hpp" 57 #include "Teuchos_StrUtils.hpp" 59 #include "Thyra_DefaultProductMultiVector.hpp" 60 #include "Thyra_DefaultProductVectorSpace.hpp" 64 using Teuchos::rcp_dynamic_cast;
71 TEUCHOS_ASSERT(blockIndex<(
int)
children_.size());
92 TEUCHOS_ASSERT(blockIndex<(
int)
children_.size());
99 TEUCHOS_ASSERT(blockIndex<(
int)
children_.size());
109 TEUCHOS_ASSERT(blockIndex<(
int)
children_.size());
117 std::stringstream ss;
119 for(
unsigned int i=0;i<
children_.size();i++) {
123 ss <<
" " <<
children_[i]->toString() <<
" ";
133 for(
unsigned int i=0;i<
children_.size();i++) {
136 int subMax =
children_[i]->LargestIndex();
137 max = max > subMax ? max : subMax;
144 Teuchos::RCP<const Thyra::LinearOpBase<double> >
146 const Teuchos::RCP<
const Thyra::BlockedLinearOpBase<double> > & blkOp)
148 return buildReorderedLinearOp(bmm,bmm,blkOp);
151 Teuchos::RCP<const Thyra::LinearOpBase<double> >
153 const Teuchos::RCP<
const Thyra::BlockedLinearOpBase<double> > & blkOp)
155 typedef RCP<const BlockReorderManager> BRMptr;
160 if(rowSz==0 && colSz==0) {
166 Teko::LinearOp linOp = blkOp->getBlock(rowLeaf.
GetIndex(),colLeaf.
GetIndex());
169 if(linOp==Teuchos::null) {
170 linOp = Thyra::zero(blkOp->productRange()->getBlock(rowLeaf.
GetIndex()),
171 blkOp->productDomain()->getBlock(colLeaf.
GetIndex()));
180 reBlkOp->beginBlockFill(1,colSz);
183 for(
int col=0;col<colSz;col++) {
184 BRMptr colPtr = colMgr.
GetBlock(col);
186 reBlkOp->setBlock(0,col,buildReorderedLinearOp(rowMgr,*colPtr,blkOp));
190 reBlkOp->endBlockFill();
198 reBlkOp->beginBlockFill(rowSz,1);
201 for(
int row=0;row<rowSz;row++) {
202 BRMptr rowPtr = rowMgr.
GetBlock(row);
204 reBlkOp->setBlock(row,0,buildReorderedLinearOp(*rowPtr,colMgr,blkOp));
208 reBlkOp->endBlockFill();
216 TEUCHOS_ASSERT(rowSz>0);
217 TEUCHOS_ASSERT(colSz>0);
220 reBlkOp->beginBlockFill(rowSz,colSz);
222 for(
int row=0;row<rowSz;row++) {
223 BRMptr rowPtr = rowMgr.
GetBlock(row);
225 for(
int col=0;col<colSz;col++) {
226 BRMptr colPtr = colMgr.
GetBlock(col);
228 reBlkOp->setBlock(row,col,buildReorderedLinearOp(*rowPtr,*colPtr,blkOp));
233 reBlkOp->endBlockFill();
260 Teuchos::RCP<const Thyra::VectorSpaceBase<double> >
262 const Teuchos::RCP<
const Thyra::ProductVectorSpaceBase<double> > & blkSpc)
264 typedef RCP<const BlockReorderManager> BRMptr;
273 return blkSpc->getBlock(leaf.
GetIndex());
276 Array<RCP<const Thyra::VectorSpaceBase<double> > > vecSpaces;
279 for(
int i=0;i<sz;i++) {
282 const RCP<const Thyra::VectorSpaceBase<double> > lvs = buildReorderedVectorSpace(*blkMgr,blkSpc);
284 vecSpaces.push_back(lvs);
288 const RCP<const Thyra::DefaultProductVectorSpace<double> > vs
289 = Thyra::productVectorSpace<double>(vecSpaces);
300 Teuchos::RCP<Thyra::MultiVectorBase<double> >
302 const Teuchos::RCP<Thyra::ProductMultiVectorBase<double> > & blkVec)
304 typedef RCP<const BlockReorderManager> BRMptr;
313 return blkVec->getNonconstMultiVectorBlock(leaf.
GetIndex());
316 Array<RCP<Thyra::MultiVectorBase<double> > > multiVecs;
317 Array<RCP<const Thyra::VectorSpaceBase<double> > > vecSpaces;
320 for(
int i=0;i<sz;i++) {
324 const RCP<const Thyra::VectorSpaceBase<double> > lvs = lmv->range();
326 multiVecs.push_back(lmv);
327 vecSpaces.push_back(lvs);
331 const RCP<const Thyra::DefaultProductVectorSpace<double> > vs
332 = Thyra::productVectorSpace<double>(vecSpaces);
335 return Thyra::defaultProductMultiVector<double>(vs,multiVecs);
343 Teuchos::RCP<const Thyra::MultiVectorBase<double> >
345 const Teuchos::RCP<
const Thyra::ProductMultiVectorBase<double> > & blkVec)
347 typedef RCP<const BlockReorderManager> BRMptr;
356 return blkVec->getMultiVectorBlock(leaf.
GetIndex());
359 Array<RCP<const Thyra::MultiVectorBase<double> > > multiVecs;
360 Array<RCP<const Thyra::VectorSpaceBase<double> > > vecSpaces;
363 for(
int i=0;i<sz;i++) {
367 const RCP<const Thyra::VectorSpaceBase<double> > lvs = lmv->range();
369 multiVecs.push_back(lmv);
370 vecSpaces.push_back(lvs);
374 const RCP<const Thyra::DefaultProductVectorSpace<double> > vs
375 = Thyra::productVectorSpace<double>(vecSpaces);
378 return Thyra::defaultProductMultiVector<double>(vs,multiVecs);
385 void buildNonconstFlatMultiVector(
const BlockReorderManager & mgr,
386 const RCP<Thyra::MultiVectorBase<double> > & blkVec,
387 Array<RCP<Thyra::MultiVectorBase<double> > > & multivecs,
388 Array<RCP<
const Thyra::VectorSpaceBase<double> > > & vecspaces)
390 int sz = mgr.GetNumBlocks();
394 const BlockReorderLeaf & leaf =
dynamic_cast<const BlockReorderLeaf &
>(mgr);
395 int index = leaf.GetIndex();
398 multivecs[index] = blkVec;
399 vecspaces[index] = blkVec->range();
402 const RCP<Thyra::ProductMultiVectorBase<double> > prodMV
403 = rcp_dynamic_cast<Thyra::ProductMultiVectorBase<double> >(blkVec);
406 for(
int i=0;i<sz;i++) {
407 const RCP<Thyra::MultiVectorBase<double> > mv = prodMV->getNonconstMultiVectorBlock(i);
408 buildNonconstFlatMultiVector(*mgr.GetBlock(i),mv,multivecs,vecspaces);
417 void buildFlatMultiVector(
const BlockReorderManager & mgr,
418 const RCP<
const Thyra::MultiVectorBase<double> > & blkVec,
419 Array<RCP<
const Thyra::MultiVectorBase<double> > > & multivecs,
420 Array<RCP<
const Thyra::VectorSpaceBase<double> > > & vecspaces)
422 int sz = mgr.GetNumBlocks();
426 const BlockReorderLeaf & leaf =
dynamic_cast<const BlockReorderLeaf &
>(mgr);
427 int index = leaf.GetIndex();
430 multivecs[index] = blkVec;
431 vecspaces[index] = blkVec->range();
434 const RCP<const Thyra::ProductMultiVectorBase<double> > prodMV
435 = rcp_dynamic_cast<
const Thyra::ProductMultiVectorBase<double> >(blkVec);
438 for(
int i=0;i<sz;i++) {
439 RCP<const Thyra::MultiVectorBase<double> > mv = prodMV->getMultiVectorBlock(i);
440 buildFlatMultiVector(*mgr.GetBlock(i),mv,multivecs,vecspaces);
450 Teuchos::RCP<Thyra::MultiVectorBase<double> >
452 const Teuchos::RCP<Thyra::ProductMultiVectorBase<double> > & blkVec)
456 Array<RCP<Thyra::MultiVectorBase<double> > > multivecs(numBlocks);
457 Array<RCP<const Thyra::VectorSpaceBase<double> > > vecspaces(numBlocks);
460 buildNonconstFlatMultiVector(mgr,blkVec,multivecs,vecspaces);
463 const RCP<Thyra::DefaultProductVectorSpace<double> > vs
464 = Thyra::productVectorSpace<double>(vecspaces);
467 return Thyra::defaultProductMultiVector<double>(vs,multivecs);
474 Teuchos::RCP<const Thyra::MultiVectorBase<double> >
476 const Teuchos::RCP<
const Thyra::ProductMultiVectorBase<double> > & blkVec)
480 Array<RCP<const Thyra::MultiVectorBase<double> > > multivecs(numBlocks);
481 Array<RCP<const Thyra::VectorSpaceBase<double> > > vecspaces(numBlocks);
484 buildFlatMultiVector(mgr,blkVec,multivecs,vecspaces);
487 const RCP<const Thyra::DefaultProductVectorSpace<double> > vs
488 = Thyra::productVectorSpace<double>(vecspaces);
491 return Thyra::defaultProductMultiVector<double>(vs,multivecs);
497 void buildFlatVectorSpace(
const BlockReorderManager & mgr,
498 const RCP<
const Thyra::VectorSpaceBase<double> > & blkSpc,
499 Array<RCP<
const Thyra::VectorSpaceBase<double> > > & vecspaces)
501 int sz = mgr.GetNumBlocks();
505 const BlockReorderLeaf & leaf =
dynamic_cast<const BlockReorderLeaf &
>(mgr);
506 int index = leaf.GetIndex();
509 vecspaces[index] = blkSpc;
512 const RCP<const Thyra::ProductVectorSpaceBase<double> > prodSpc
513 = rcp_dynamic_cast<
const Thyra::ProductVectorSpaceBase<double> >(blkSpc);
516 for(
int i=0;i<sz;i++) {
517 RCP<const Thyra::VectorSpaceBase<double> > space = prodSpc->getBlock(i);
518 buildFlatVectorSpace(*mgr.GetBlock(i),space,vecspaces);
526 Teuchos::RCP<const Thyra::VectorSpaceBase<double> >
527 buildFlatVectorSpace(
const BlockReorderManager & mgr,
528 const Teuchos::RCP<
const Thyra::VectorSpaceBase<double> > & blkSpc)
530 int numBlocks = mgr.LargestIndex()+1;
532 Array<RCP<const Thyra::VectorSpaceBase<double> > > vecspaces(numBlocks);
535 buildFlatVectorSpace(mgr,blkSpc,vecspaces);
538 return Thyra::productVectorSpace<double>(vecspaces);
548 static void tokenize(std::string srcInput,std::string whitespace,std::string prefer,
549 std::vector<std::string> & tokens)
551 std::string input = srcInput;
552 std::vector<std::string> wsTokens;
553 std::size_t endPos = input.length()-1;
554 while(endPos<input.length()) {
555 std::size_t next = input.find_first_of(whitespace);
560 if(next!=std::string::npos) {
561 s = input.substr(0,next);
564 input = input.substr(next+1,endPos);
571 endPos = input.length()-1;
575 wsTokens.push_back(s);
578 for(
unsigned int i=0;i<wsTokens.size();i++) {
582 std::size_t endPos = input.length()-1;
583 while(endPos<input.length()) {
584 std::size_t next = input.find_first_of(prefer);
586 std::string s = input;
587 if(next>0 && next<input.length()) {
590 s = input.substr(0,next);
592 input = input.substr(next,endPos);
596 s = input.substr(0,next+1);
598 input = input.substr(next+1,endPos);
603 endPos = input.length()-1;
613 static std::vector<std::string>::const_iterator
614 buildSubBlock(std::vector<std::string>::const_iterator begin,
615 std::vector<std::string>::const_iterator end, std::vector<std::string> & subBlock)
617 std::stack<std::string> matched;
618 std::vector<std::string>::const_iterator itr;
619 for(itr=begin;itr!=end;++itr) {
621 subBlock.push_back(*itr);
624 if(*itr==
"[") matched.push(
"[");
625 else if(*itr==
"]") matched.pop();
632 TEUCHOS_ASSERT(matched.empty());
638 static RCP<BlockReorderManager> blockedReorderFromTokens(
const std::vector<std::string> & tokens)
645 TEUCHOS_ASSERT(*(tokens.begin())==
"[")
646 TEUCHOS_ASSERT(*(tokens.end()-1)=="]");
648 std::vector<RCP<
Teko::BlockReorderManager> > vecRMgr;
649 std::vector<std::
string>::const_iterator itr = tokens.begin()+1;
650 while(itr!=tokens.end()-1) {
652 std::vector<std::string> subBlock;
653 itr = buildSubBlock(itr,tokens.end()-1,subBlock);
656 vecRMgr.push_back(blockedReorderFromTokens(subBlock));
664 for(
unsigned int i=0;i<vecRMgr.size();i++)
665 rMgr->SetBlock(i,vecRMgr[i]);
686 std::vector<std::string> tokens;
691 tokenize(reorder,
" \t\n",
"[]",tokens);
694 Teuchos::RCP<BlockReorderManager> mgr = blockedReorderFromTokens(tokens);
Teuchos::RCP< Thyra::MultiVectorBase< double > > buildReorderedMultiVector(const BlockReorderManager &mgr, const Teuchos::RCP< Thyra::ProductMultiVectorBase< double > > &blkVec)
Convert a flat multi vector into a reordered multivector.
std::vector< Teuchos::RCP< BlockReorderManager > > children_
Definitions of the subblocks.
Class that describes how a flat blocked operator should be reordered.
virtual const Teuchos::RCP< BlockReorderManager > GetBlock(int blockIndex)
Get a particular block. If there is no block at this index location return a new one.
Teuchos::RCP< const Thyra::VectorSpaceBase< double > > buildReorderedVectorSpace(const BlockReorderManager &mgr, const Teuchos::RCP< const Thyra::ProductVectorSpaceBase< double > > &blkSpc)
Use the BlockReorderManager to change a flat vector space into a composite vector space...
virtual int LargestIndex() const
Largest index in this manager.
int GetIndex() const
Get the the index that is stored in this block.
BlockReorderManager()
Basic empty constructor.
BlockedLinearOp createBlockedOp()
Build a new blocked linear operator.
virtual std::string toString() const
For sanities sake, print a readable string.
virtual int GetNumBlocks() const
Gets the number of subblocks.
virtual void SetBlock(int blockIndex, int reorder)
Sets the sublock to a specific index value.
Teuchos::RCP< const BlockReorderManager > blockedReorderFromString(std::string &reorder)
Convert a string to a block reorder manager object.