42 #include "Thyra_EpetraOperatorWrapper.hpp" 44 #include "Thyra_DetachedSpmdVectorView.hpp" 45 #include "Thyra_DefaultProductVector.hpp" 46 #include "Thyra_ProductVectorSpaceBase.hpp" 47 #include "Thyra_SpmdVectorBase.hpp" 48 #include "Thyra_EpetraLinearOp.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" 71 : useTranspose_(false),
73 range_(thyraOp->range()),
74 domain_(thyraOp->domain()),
75 comm_(getEpetraComm(*thyraOp)),
78 label_(thyraOp->description())
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 =
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;
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 =
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> >
164 opRange = ( !useTranspose_ ? range_ : domain_ ),
165 opDomain = ( !useTranspose_ ? domain_ : range_ );
167 const RCP<VectorBase<double> > tx = createMember(opDomain);
170 const RCP<VectorBase<double> > ty = createMember(opRange);
172 Thyra::apply<double>( *thyraOp_, !useTranspose_ ?
NOTRANS :
CONJTRANS, *tx, ty.ptr());
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 =
217 if (nonnull(prod_vs))
220 const RCP<const SpmdVectorSpaceBase<double> > spmd_vs =
231 Teuchos::RCP<const Thyra::LinearOpBase<double> >
232 Thyra::makeEpetraWrapper(
const RCP<
const LinearOpBase<double> > &thyraOp)
235 Teuchos::rcp(
new EpetraOperatorWrapper(thyraOp))
Create an explicit detached mutable (non-const) view of all of the local elements on this process of ...
virtual RCP< const VectorSpaceBase< Scalar > > range() const =0
Return a smart pointer for the range space for this operator.
EpetraOperatorWrapper(const RCP< const LinearOpBase< double > > &thyraOp)
virtual Teuchos::RCP< const VectorSpaceBase< Scalar > > getBlock(const int k) const =0
Returns a vector space for the kth (zero-based) block.
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...
Use the non-transposed operator.
int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Use the transposed operator with complex-conjugate clements (same as TRANS for real scalar types)...
void copyThyraIntoEpetra(const VectorBase< double > &thyraVec, Epetra_MultiVector &x) const
Abstract interface for finite-dimensional dense vectors.
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 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
Create an explicit detached non-mutable (const) view of all of the local elements on this process of ...
int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Base abstract VectorSpaceBase class for all SPMD-based vector spaces.