44 #include "Thyra_DetachedSpmdVectorView.hpp" 45 #include "Thyra_DefaultProductVector.hpp" 46 #include "Thyra_ProductVectorSpaceBase.hpp" 47 #include "Thyra_SpmdVectorBase.hpp" 51 # include "Epetra_MpiComm.h" 53 #include "Epetra_SerialComm.h" 54 #include "Epetra_Vector.h" 57 # include "Teuchos_DefaultMpiComm.hpp" 59 #include "Teuchos_DefaultSerialComm.hpp" 69 const RCP<
const LinearOpBase<double> > &thyraOp
71 : useTranspose_(false),
73 range_(thyraOp->range()),
74 domain_(thyraOp->domain()),
75 comm_(getEpetraComm(*thyraOp)),
78 label_(thyraOp->description())
83 const Ptr<VectorBase<double> > &thyraVec)
const 86 using Teuchos::rcpFromPtr;
87 using Teuchos::rcp_dynamic_cast;
89 const int numVecs = x.NumVectors();
91 TEUCHOS_TEST_FOR_EXCEPTION(numVecs != 1, std::runtime_error,
92 "epetraToThyra does not work with MV dimension != 1");
94 const RCP<ProductVectorBase<double> > prodThyraVec =
95 castOrCreateNonconstProductVectorBase(rcpFromPtr(thyraVec));
97 const ArrayView<const double> epetraData(x[0], x.Map().NumMyElements());
103 const int numBlocks = prodThyraVec->productSpace()->numBlocks();
104 for (
int b = 0; b < numBlocks; ++b) {
105 const RCP<VectorBase<double> > vec_b = prodThyraVec->getNonconstVectorBlock(b);
106 const RCP<const SpmdVectorSpaceBase<double> > spmd_vs_b =
107 rcp_dynamic_cast<
const SpmdVectorSpaceBase<double> >(vec_b->space(),
true);
108 DetachedSpmdVectorView<double> view(vec_b);
109 const ArrayRCP<double> thyraData = view.sv().values();
110 const int localNumElems = spmd_vs_b->localSubDim();
111 for (
int i=0; i < localNumElems; ++i) {
112 thyraData[i] = epetraData[i+offset];
114 offset += localNumElems;
121 const VectorBase<double>& thyraVec, Epetra_MultiVector& x)
const 124 using Teuchos::rcpFromRef;
125 using Teuchos::rcp_dynamic_cast;
127 const int numVecs = x.NumVectors();
129 TEUCHOS_TEST_FOR_EXCEPTION(numVecs != 1, std::runtime_error,
130 "epetraToThyra does not work with MV dimension != 1");
132 const RCP<const ProductVectorBase<double> > prodThyraVec =
133 castOrCreateProductVectorBase(rcpFromRef(thyraVec));
135 const ArrayView<double> epetraData(x[0], x.Map().NumMyElements());
139 const int numBlocks = prodThyraVec->productSpace()->numBlocks();
140 for (
int b = 0; b < numBlocks; ++b) {
141 const RCP<const VectorBase<double> > vec_b = prodThyraVec->getVectorBlock(b);
142 const RCP<const SpmdVectorSpaceBase<double> > spmd_vs_b =
143 rcp_dynamic_cast<
const SpmdVectorSpaceBase<double> >(vec_b->space(),
true);
144 ConstDetachedSpmdVectorView<double> view(vec_b);
145 const ArrayRCP<const double> thyraData = view.sv().values();
146 const int localNumElems = spmd_vs_b->localSubDim();
147 for (
int i=0; i < localNumElems; ++i) {
148 epetraData[i+offset] = thyraData[i];
150 offset += localNumElems;
160 Epetra_MultiVector& Y)
const 163 const RCP<const VectorSpaceBase<double> >
167 const RCP<VectorBase<double> > tx = createMember(opDomain);
170 const RCP<VectorBase<double> > ty = createMember(opRange);
182 Epetra_MultiVector& Y)
const 184 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::runtime_error,
185 "EpetraOperatorWrapper::ApplyInverse not implemented");
192 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::runtime_error,
193 "EpetraOperatorWrapper::NormInf not implemated");
201 RCP<const Epetra_Comm>
206 using Teuchos::rcp_dynamic_cast;
207 using Teuchos::SerialComm;
209 using Teuchos::MpiComm;
212 RCP<const VectorSpaceBase<double> > vs = thyraOp.range();
214 const RCP<const ProductVectorSpaceBase<double> > prod_vs =
215 rcp_dynamic_cast<
const ProductVectorSpaceBase<double> >(vs);
217 if (nonnull(prod_vs))
218 vs = prod_vs->getBlock(0);
220 const RCP<const SpmdVectorSpaceBase<double> > spmd_vs =
221 rcp_dynamic_cast<
const SpmdVectorSpaceBase<double> >(vs,
true);
231 Teuchos::RCP<const Thyra::LinearOpBase<double> >
232 Thyra::makeEpetraWrapper(
const RCP<
const LinearOpBase<double> > &thyraOp)
235 Teuchos::rcp(
new EpetraOperatorWrapper(thyraOp))
EpetraOperatorWrapper(const RCP< const LinearOpBase< double > > &thyraOp)
RCP< const Epetra_Map > get_Epetra_Map(const VectorSpaceBase< double > &vs, const RCP< const Epetra_Comm > &comm)
Get (or create) an Epetra_Map object given an VectorSpaceBase object an optionally an extra Epetra_Co...
int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
static RCP< const Epetra_Comm > getEpetraComm(const LinearOpBase< double > &thyraOp)
void copyThyraIntoEpetra(const VectorBase< double > &thyraVec, Epetra_MultiVector &x) const
RCP< const EpetraLinearOp > epetraLinearOp(const RCP< const Epetra_Operator > &op, EOpTransp opTrans=NOTRANS, EApplyEpetraOpAs applyAs=EPETRA_OP_APPLY_APPLY, EAdjointEpetraOp adjointSupport=EPETRA_OP_ADJOINT_SUPPORTED, const RCP< const VectorSpaceBase< double > > &range=Teuchos::null, const RCP< const VectorSpaceBase< double > > &domain=Teuchos::null)
Dynamically allocate a nonconst EpetraLinearOp to wrap a const Epetra_Operator object.
RCP< const VectorSpaceBase< double > > range_
RCP< const VectorSpaceBase< double > > domain_
RCP< const Epetra_Comm > get_Epetra_Comm(const Teuchos::Comm< Ordinal > &comm)
Get (or create) and Epetra_Comm given a Teuchos::Comm object.
void copyEpetraIntoThyra(const Epetra_MultiVector &x, const Ptr< VectorBase< double > > &thyraVec) const
RCP< const LinearOpBase< double > > thyraOp_
int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const