47 #include "Teko_TpetraStridedMappingStrategy.hpp"
48 #include "Teko_InterlacedTpetra.hpp"
49 #include "Teko_TpetraHelpers.hpp"
51 #include "Thyra_TpetraThyraWrappers.hpp"
52 #include "Thyra_TpetraLinearOp.hpp"
53 #include "Thyra_DefaultProductMultiVector.hpp"
54 #include "Thyra_DefaultProductVectorSpace.hpp"
55 #include "Thyra_DefaultSpmdMultiVector.hpp"
56 #include "Thyra_DefaultBlockedLinearOp.hpp"
60 using Teuchos::rcp_dynamic_cast;
63 namespace TpetraHelpers {
75 TpetraStridedMappingStrategy::TpetraStridedMappingStrategy(
76 const std::vector<int>& vars,
const RCP<
const Tpetra::Map<LO, GO, NT> >& map,
77 const Teuchos::Comm<int>& comm) {
80 buildBlockTransferData(vars, rangeMap_, comm);
91 void TpetraStridedMappingStrategy::copyTpetraIntoThyra(
92 const Tpetra::MultiVector<ST, LO, GO, NT>& X,
93 const Teuchos::Ptr<Thyra::MultiVectorBase<ST> >& thyra_X)
const {
94 int count = X.getNumVectors();
96 std::vector<RCP<Tpetra::MultiVector<ST, LO, GO, NT> > > subX;
99 Strided::buildSubVectors(blockMaps_, subX, count);
102 Strided::one2many(subX, X, blockImport_);
105 Teuchos::Ptr<Thyra::ProductMultiVectorBase<ST> > prod_X =
106 Teuchos::ptr_dynamic_cast<Thyra::ProductMultiVectorBase<ST> >(thyra_X);
107 for (
unsigned int i = 0; i < blockMaps_.size(); i++) {
108 RCP<Thyra::TpetraMultiVector<ST, LO, GO, NT> > vec =
109 rcp_dynamic_cast<Thyra::TpetraMultiVector<ST, LO, GO, NT> >(
110 prod_X->getNonconstMultiVectorBlock(i),
true);
112 fillDefaultSpmdMultiVector(vec, subX[i]);
124 void TpetraStridedMappingStrategy::copyThyraIntoTpetra(
125 const RCP<
const Thyra::MultiVectorBase<ST> >& thyra_Y,
126 Tpetra::MultiVector<ST, LO, GO, NT>& Y)
const {
127 std::vector<RCP<const Tpetra::MultiVector<ST, LO, GO, NT> > > subY;
128 RCP<const Thyra::DefaultProductMultiVector<ST> > prod_Y =
129 rcp_dynamic_cast<
const Thyra::DefaultProductMultiVector<ST> >(thyra_Y);
132 for (
unsigned int i = 0; i < blockMaps_.size(); i++) {
133 RCP<const Thyra::TpetraMultiVector<ST, LO, GO, NT> > tmv =
134 rcp_dynamic_cast<
const Thyra::TpetraMultiVector<ST, LO, GO, NT> >(
135 prod_Y->getMultiVectorBlock(i),
true);
136 subY.push_back(tmv->getConstTpetraMultiVector());
140 Strided::associateSubVectors(blockMaps_, subY);
143 Strided::many2one(Y, subY, blockExport_);
158 void TpetraStridedMappingStrategy::buildBlockTransferData(
159 const std::vector<int>& vars,
const Teuchos::RCP<
const Tpetra::Map<LO, GO, NT> >& baseMap,
160 const Teuchos::Comm<int>& comm) {
162 Strided::buildSubMaps(*baseMap, vars, comm, blockMaps_);
163 Strided::buildExportImport(*baseMap, blockMaps_, blockExport_, blockImport_);
176 const Teuchos::RCP<Thyra::BlockedLinearOpBase<ST> >
177 TpetraStridedMappingStrategy::buildBlockedThyraOp(
178 const RCP<
const Tpetra::CrsMatrix<ST, LO, GO, NT> >& crsContent,
179 const std::string& label)
const {
180 int dim = blockMaps_.size();
182 RCP<Thyra::DefaultBlockedLinearOp<ST> > A = Thyra::defaultBlockedLinearOp<ST>();
184 A->beginBlockFill(dim, dim);
185 for (
int i = 0; i < dim; i++) {
186 for (
int j = 0; j < dim; j++) {
188 std::stringstream ss;
189 ss << label <<
"_" << i <<
"," << j;
192 RCP<Tpetra::CrsMatrix<ST, LO, GO, NT> > Aij =
193 Strided::buildSubBlock(i, j, crsContent, blockMaps_);
194 A->setNonconstBlock(i, j,
195 Thyra::tpetraLinearOp<ST, LO, GO, NT>(
196 Thyra::tpetraVectorSpace<ST, LO, GO, NT>(Aij->getRangeMap()),
197 Thyra::tpetraVectorSpace<ST, LO, GO, NT>(Aij->getDomainMap()), Aij));
215 void TpetraStridedMappingStrategy::rebuildBlockedThyraOp(
216 const RCP<
const Tpetra::CrsMatrix<ST, LO, GO, NT> >& crsContent,
217 const RCP<Thyra::BlockedLinearOpBase<ST> >& A)
const {
218 int dim = blockMaps_.size();
220 for (
int i = 0; i < dim; i++) {
221 for (
int j = 0; j < dim; j++) {
223 RCP<Thyra::LinearOpBase<ST> > Aij = A->getNonconstBlock(i, j);
224 RCP<Thyra::TpetraLinearOp<ST, LO, GO, NT> > tAij =
225 rcp_dynamic_cast<Thyra::TpetraLinearOp<ST, LO, GO, NT> >(Aij,
true);
226 RCP<Tpetra::CrsMatrix<ST, LO, GO, NT> > eAij =
227 rcp_dynamic_cast<Tpetra::CrsMatrix<ST, LO, GO, NT> >(tAij->getTpetraOperator(),
true);
230 Strided::rebuildSubBlock(i, j, crsContent, blockMaps_, *eAij);