47 #include "Teko_TpetraThyraConverter.hpp"
48 #include "Tpetra_Core.hpp"
51 #include "Teuchos_Array.hpp"
52 #include "Teuchos_ArrayRCP.hpp"
55 #include "Thyra_DefaultProductVectorSpace.hpp"
56 #include "Thyra_DefaultProductMultiVector.hpp"
57 #include "Thyra_SpmdMultiVectorBase.hpp"
58 #include "Thyra_SpmdVectorSpaceBase.hpp"
59 #include "Thyra_MultiVectorStdOps.hpp"
66 using Teuchos::ptr_dynamic_cast;
69 using Teuchos::rcp_dynamic_cast;
70 using Teuchos::rcpFromRef;
73 namespace TpetraHelpers {
79 void blockTpetraToThyra(
int numVectors, Teuchos::ArrayRCP<const ST> tpetraData,
int leadingDim,
80 const Teuchos::Ptr<Thyra::MultiVectorBase<ST> >& mv,
int& localDim) {
84 const Ptr<Thyra::ProductMultiVectorBase<ST> > prodMV =
85 ptr_dynamic_cast<Thyra::ProductMultiVectorBase<ST> >(mv);
86 if (prodMV == Teuchos::null) {
88 const Ptr<Thyra::SpmdMultiVectorBase<ST> > spmdX =
89 ptr_dynamic_cast<Thyra::SpmdMultiVectorBase<ST> >(mv,
true);
90 const RCP<const Thyra::SpmdVectorSpaceBase<ST> > spmdVS = spmdX->spmdSpace();
92 int localSubDim = spmdVS->localSubDim();
94 Thyra::Ordinal thyraLeadingDim = 0;
96 Teuchos::ArrayRCP<ST> thyraData_arcp;
97 Teuchos::ArrayView<ST> thyraData;
98 spmdX->getNonconstLocalData(Teuchos::outArg(thyraData_arcp), Teuchos::outArg(thyraLeadingDim));
99 thyraData = thyraData_arcp();
101 for (
int i = 0; i < localSubDim; i++) {
103 for (
int v = 0; v < numVectors; v++)
104 thyraData[i + thyraLeadingDim * v] = tpetraData[i + leadingDim * v];
108 localDim = localSubDim;
114 Teuchos::ArrayRCP<const ST> localData = tpetraData;
117 for (
int blkIndex = 0; blkIndex < prodMV->productSpace()->numBlocks(); blkIndex++) {
119 const RCP<Thyra::MultiVectorBase<ST> > blockVec = prodMV->getNonconstMultiVectorBlock(blkIndex);
122 blockTpetraToThyra(numVectors, localData, leadingDim, blockVec.ptr(), subDim);
136 void blockTpetraToThyra(
const Tpetra::MultiVector<ST, LO, GO, NT>& tpetraX,
137 const Teuchos::Ptr<Thyra::MultiVectorBase<ST> >& thyraX) {
138 TEUCHOS_ASSERT((Tpetra::global_size_t)thyraX->range()->dim() == tpetraX.getGlobalLength());
141 LO leadingDim = 0, localDim = 0;
142 leadingDim = tpetraX.getStride();
143 Teuchos::ArrayRCP<const ST> tpetraData = tpetraX.get1dView();
145 int numVectors = tpetraX.getNumVectors();
147 blockTpetraToThyra(numVectors, tpetraData, leadingDim, thyraX.ptr(), localDim);
149 TEUCHOS_ASSERT((
size_t)localDim == tpetraX.getLocalLength());
152 void blockThyraToTpetra(LO numVectors, Teuchos::ArrayRCP<ST> tpetraData, LO leadingDim,
153 const Teuchos::RCP<
const Thyra::MultiVectorBase<ST> >& tX, LO& localDim) {
157 const RCP<const Thyra::ProductMultiVectorBase<ST> > prodX =
158 rcp_dynamic_cast<
const Thyra::ProductMultiVectorBase<ST> >(tX);
159 if (prodX == Teuchos::null) {
163 RCP<const Thyra::SpmdMultiVectorBase<ST> > spmdX =
164 rcp_dynamic_cast<
const Thyra::SpmdMultiVectorBase<ST> >(tX,
true);
165 RCP<const Thyra::SpmdVectorSpaceBase<ST> > spmdVS = spmdX->spmdSpace();
167 Thyra::Ordinal thyraLeadingDim = 0;
168 Teuchos::ArrayView<const ST> thyraData;
169 Teuchos::ArrayRCP<const ST> thyraData_arcp;
170 spmdX->getLocalData(Teuchos::outArg(thyraData_arcp), Teuchos::outArg(thyraLeadingDim));
171 thyraData = thyraData_arcp();
173 LO localSubDim = spmdVS->localSubDim();
174 for (LO i = 0; i < localSubDim; i++) {
176 for (LO v = 0; v < numVectors; v++) {
177 tpetraData[i + leadingDim * v] = thyraData[i + thyraLeadingDim * v];
182 localDim = localSubDim;
187 const RCP<const Thyra::ProductVectorSpaceBase<ST> > prodVS = prodX->productSpace();
190 Teuchos::ArrayRCP<ST> localData = tpetraData;
193 for (
int blkIndex = 0; blkIndex < prodVS->numBlocks(); blkIndex++) {
197 blockThyraToTpetra(numVectors, localData, leadingDim, prodX->getMultiVectorBlock(blkIndex),
215 void blockThyraToTpetra(
const Teuchos::RCP<
const Thyra::MultiVectorBase<ST> >& thyraX,
216 Tpetra::MultiVector<ST, LO, GO, NT>& tpetraX) {
218 LO numVectors = thyraX->domain()->dim();
221 TEUCHOS_ASSERT((
size_t)numVectors == tpetraX.getNumVectors());
222 TEUCHOS_ASSERT((Tpetra::global_size_t)thyraX->range()->dim() == tpetraX.getGlobalLength());
225 LO leadingDim = 0, localDim = 0;
226 leadingDim = tpetraX.getStride();
227 Teuchos::ArrayRCP<ST> tpetraData = tpetraX.get1dViewNonConst();
230 blockThyraToTpetra(numVectors, tpetraData, leadingDim, thyraX, localDim);
233 TEUCHOS_ASSERT((
size_t)localDim == tpetraX.getLocalLength());
236 void thyraVSToTpetraMap(std::vector<GO>& myIndicies,
int blockOffset,
237 const Thyra::VectorSpaceBase<ST>& vs,
int& localDim) {
241 const RCP<const Thyra::ProductVectorSpaceBase<ST> > prodVS =
242 rcp_dynamic_cast<
const Thyra::ProductVectorSpaceBase<ST> >(rcpFromRef(vs));
245 if (prodVS == Teuchos::null) {
249 const RCP<const Thyra::SpmdVectorSpaceBase<ST> > spmdVS =
250 rcp_dynamic_cast<
const Thyra::SpmdVectorSpaceBase<ST> >(rcpFromRef(vs));
251 TEUCHOS_TEST_FOR_EXCEPTION(spmdVS == Teuchos::null, std::runtime_error,
252 "thyraVSToTpetraMap requires all subblocks to be SPMD");
255 int localOffset = spmdVS->localOffset();
256 int localSubDim = spmdVS->localSubDim();
259 for (
int i = 0; i < localSubDim; i++) myIndicies.push_back(blockOffset + localOffset + i);
261 localDim += localSubDim;
267 for (
int blkIndex = 0; blkIndex < prodVS->numBlocks(); blkIndex++) {
271 thyraVSToTpetraMap(myIndicies, blockOffset, *prodVS->getBlock(blkIndex), subDim);
273 blockOffset += prodVS->getBlock(blkIndex)->dim();
281 const RCP<Tpetra::Map<LO, GO, NT> > thyraVSToTpetraMap(
282 const Thyra::VectorSpaceBase<ST>& vs,
283 const RCP<
const Teuchos::Comm<Thyra::Ordinal> >& ) {
285 std::vector<GO> myGIDs;
288 thyraVSToTpetraMap(myGIDs, 0, vs, localDim);
290 TEUCHOS_ASSERT(myGIDs.size() == (size_t)localDim);
296 return rcp(
new Tpetra::Map<LO, GO, NT>(vs.dim(), Teuchos::ArrayView<const GO>(myGIDs), 0,
297 Tpetra::getDefaultComm()));