43 #include "Epetra_Vector.h" 44 #include "Epetra_Map.h" 55 ,
const Teuchos::RCP<const Epetra_Operator> op[]
56 ,
const Teuchos::ETransp opTrans[]
65 ,
const Teuchos::RCP<const Epetra_Operator> op[]
66 ,
const Teuchos::ETransp opTrans[]
71 TEUCHOS_TEST_FOR_EXCEPTION(
72 numOp < 1, std::invalid_argument
73 ,
"ProductOperator::initialize(...): Error!" 80 std::copy( op, op + numOp,
Op_.begin() );
81 std::copy( opTrans, opTrans + numOp,
Op_trans_.begin() );
82 std::copy( opInverse, opInverse + numOp,
Op_inverse_.begin() );
91 ,Teuchos::RCP<const Epetra_Operator> op[]
92 ,Teuchos::ETransp opTrans[]
97 TEUCHOS_TEST_FOR_EXCEPTION(
98 (op != NULL || opTrans != NULL || opInverse!=NULL) && numOp==NULL
99 ,std::invalid_argument
100 ,
"ProductOperator::uninitialize(...): Error!" 105 if(op) std::copy(
Op_.begin(),
Op_.end(), op );
119 ,Teuchos::ETransp opTrans
121 ,
const Epetra_MultiVector &X_k
122 ,Epetra_MultiVector *Y_k
125 Epetra_Operator &Op_k =
const_cast<Epetra_Operator&
>(*
Op_[k]);
126 bool oldUseTranspose = Op_k.UseTranspose();
127 Op_k.SetUseTranspose((opTrans==Teuchos::NO_TRANS) != (
Op_trans_[k]==Teuchos::NO_TRANS));
129 const int err = ! applyInverse_k ?
Op_[k]->Apply(X_k,*Y_k) :
Op_[k]->ApplyInverse(X_k,*Y_k);
130 Op_k.SetUseTranspose(oldUseTranspose);
131 TEUCHOS_TEST_FOR_EXCEPTION(
132 err!=0, std::runtime_error,
"ProductOperator::applyConstituent(...): Error," 133 " Op["<<k<<
"]." << (!applyInverse_k?
"Apply":
"ApplyInverse") <<
"(...) " 134 "returned err = " << err <<
" with Op["<<k<<
"].UseTranspose() = "<<
150 const int numOp = this->
num_Op();
158 for(
int k = numOp-1; k >= 0; --k ) {
159 const Epetra_MultiVector &X_k = ( k==numOp-1 ? X : *
range_vecs_[k] );
160 Epetra_MultiVector &Y_k = ( k==0 ? Y : *
range_vecs_[k-1] );
168 for(
int k = 0; k <= numOp-1; ++k ) {
169 const Epetra_MultiVector &X_k = ( k==0 ? X : *
domain_vecs_[k-1] );
170 Epetra_MultiVector &Y_k = ( k==numOp-1 ? Y : *
domain_vecs_[k] );
180 const int numOp = this->
num_Op();
188 for(
int k = 0; k <= numOp-1; ++k ) {
189 const Epetra_MultiVector &X_k = ( k==0 ? X : *
domain_vecs_[k-1] );
190 Epetra_MultiVector &Y_k = ( k==numOp-1 ? Y : *
domain_vecs_[k] );
198 for(
int k = numOp-1; k >= 0; --k ) {
199 const Epetra_MultiVector &X_k = ( k==numOp-1 ? X : *
range_vecs_[k] );
200 Epetra_MultiVector &Y_k = ( k==0 ? Y : *
range_vecs_[k-1] );
235 return Op_.front()->OperatorRangeMap().Comm();
242 return (
Op_trans_.back()==Teuchos::NO_TRANS
243 ?
Op_.back()->OperatorDomainMap()
244 :
Op_.back()->OperatorRangeMap()
252 return (
Op_trans_.front()==Teuchos::NO_TRANS
253 ?
Op_.front()->OperatorRangeMap()
254 :
Op_.front()->OperatorDomainMap()
262 const int numOp = this->
num_Op ();
283 for (
int k = numOp-1; k >= 1; --k) {
308 for (
int k = 0; k <= numOp - 2; ++k) {
void applyConstituent(const int k, Teuchos::ETransp Op_trans, EApplyMode Op_inverse, const Epetra_MultiVector &X_k, Epetra_MultiVector *Y_k) const
Apply the kth aggregate operator M[k] correctly.
const Epetra_Map & OperatorDomainMap() const
const Epetra_Comm & Comm() const
int SetUseTranspose(bool UseTranspose)
void assertInitialized() const
void initialize(const int num_Op, const Teuchos::RCP< const Epetra_Operator > Op[], const Teuchos::ETransp Op_trans[], const EApplyMode Op_inverse[])
Setup with constituent operators.
ProductOperator()
Construct to uninitialized.
EpetraExt::BlockCrsMatrix: A class for constructing a distributed block matrix.
int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
int num_Op() const
Return the number of aggregate opeators.
int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
void initializeTempVecs(bool applyInverse) const
void uninitialize(int *num_Op, Teuchos::RCP< const Epetra_Operator > Op[], Teuchos::ETransp Op_trans[], EApplyMode p_inverse[])
Set to an uninitialized state and wipe out memory.
const char * Label() const
const Epetra_Map & OperatorRangeMap() const
bool UseTranspose() const