Xpetra Version of the Day
Loading...
Searching...
No Matches
Xpetra_EpetraOperator.hpp
Go to the documentation of this file.
1// @HEADER
2//
3// ***********************************************************************
4//
5// Xpetra: A linear algebra interface package
6// Copyright 2012 Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact
39// Jonathan Hu (jhu@sandia.gov)
40// Andrey Prokopenko (aprokop@sandia.gov)
41// Ray Tuminaro (rstumin@sandia.gov)
42//
43// ***********************************************************************
44//
45// @HEADER
46#ifndef XPETRA_EPETRAOPERATOR_HPP
47#define XPETRA_EPETRAOPERATOR_HPP
48
50
51#include <Epetra_Operator.h>
52#include <Epetra_Map.h>
53
54#include "Xpetra_Map.hpp"
55#include "Xpetra_EpetraMap.hpp"
56#include "Xpetra_MultiVector.hpp"
58#include "Xpetra_Operator.hpp"
59
60#include "Xpetra_Utils.hpp"
61
62namespace Xpetra {
63
64 template<class EpetraGlobalOrdinal, class Node>
65 class EpetraOperator : public Operator<double, int, EpetraGlobalOrdinal,Node>
66 {
67 typedef double Scalar;
68 typedef int LocalOrdinal;
69 typedef EpetraGlobalOrdinal GlobalOrdinal;
70
71 public:
73
74
77 {
78 XPETRA_MONITOR("EpetraOperator::getDomainMap()");
79 return toXpetra<GlobalOrdinal,Node>(op_->OperatorDomainMap());
80 }
81
82
85 {
86 XPETRA_MONITOR("EpetraOperator::getRangeMap()");
87 return toXpetra<GlobalOrdinal,Node>(op_->OperatorRangeMap());
88 }
89
90
92
97 virtual void
103 {
104 XPETRA_MONITOR("EpetraOperator::apply");
105
106 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, X, eX, "Xpetra::EpetraOperator->apply(): cannot cast input to Xpetra::EpetraMultiVectorT");
107 XPETRA_DYNAMIC_CAST( EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, Y, eY, "Xpetra::EpetraOperator->apply(): cannot cast input to Xpetra::EpetraMultiVectorT");
108
110 "Xpetra::EpetraOperator->apply(): can only accept mode == NO_TRANS or mode == TRANS");
112 "Xpetra::EpetraOperator->apply(): cannot apply transpose as underlying Epetra operator does not support it");
113
114 // Helper vector for string A*X
115 RCP<Epetra_MultiVector> epY = eY.getEpetra_MultiVector();
117 tmp->PutScalar(0.0);
118
119 op_->SetUseTranspose(mode == Teuchos::TRANS);
120 XPETRA_ERR_CHECK(op_->Apply(*eX.getEpetra_MultiVector(), *tmp));
121
122 // calculate alpha * A * x + beta * y
123 XPETRA_ERR_CHECK(eY.getEpetra_MultiVector()->Update(alpha, *tmp, beta));
124 }
125
126
128 virtual bool hasTransposeApply() const
129 {
130 // We do not currently use transpose, try setting it
131 int err = op_->SetUseTranspose(true);
132 op_->SetUseTranspose(false);
133 return (err == 0);
134 }
135
137
139
140
142 std::string description() const { XPETRA_MONITOR("EpetraOperator::description"); return "Epetra_Operator"; }
143
146 {
147 XPETRA_MONITOR("EpetraOperator::describe");
148 out << "Epetra_Operator" << std::endl;
149 }
150
152
154
155
157 EpetraOperator(const Teuchos::RCP<Epetra_Operator> &op) : op_(op) { } //TODO removed const
158
167
169
170 private:
173
174 }; // EpetraOperator class
175
176
177 template<class EpetraGlobalOrdinal, class Node>
178 class EpetraInverseOperator : public Operator<double, int, EpetraGlobalOrdinal,Node>
179 {
180 typedef double Scalar;
181 typedef int LocalOrdinal;
182 typedef EpetraGlobalOrdinal GlobalOrdinal;
183
184 public:
186
189 XPETRA_MONITOR("EpetraOperator::getDomainMap()");
190 return toXpetra<GlobalOrdinal,Node>(op_->OperatorDomainMap());
191 }
192
195 XPETRA_MONITOR("EpetraOperator::getRangeMap()");
196 return toXpetra<GlobalOrdinal,Node>(op_->OperatorRangeMap());
197 }
198
200
205 virtual void
211 XPETRA_MONITOR("EpetraOperator::apply");
212
213 XPETRA_DYNAMIC_CAST(const EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, X, eX, "Xpetra::EpetraOperator->apply(): cannot cast input to Xpetra::EpetraMultiVectorT");
214 XPETRA_DYNAMIC_CAST( EpetraMultiVectorT<GlobalOrdinal XPETRA_COMMA Node>, Y, eY, "Xpetra::EpetraOperator->apply(): cannot cast input to Xpetra::EpetraMultiVectorT");
215
217 "Xpetra::EpetraOperator->apply(): can only accept mode == NO_TRANS or mode == TRANS");
219 "Xpetra::EpetraOperator->apply(): cannot apply transpose as underlying Epetra operator does not support it");
220
221 // Helper vector for string A*X
222 RCP<Epetra_MultiVector> epY = eY.getEpetra_MultiVector();
224 tmp->PutScalar(0.0);
225
226 op_->SetUseTranspose(mode == Teuchos::TRANS);
227 XPETRA_ERR_CHECK(op_->ApplyInverse(*eX.getEpetra_MultiVector(), *tmp));
228
229 // calculate alpha * A * x + beta * y
230 XPETRA_ERR_CHECK(eY.getEpetra_MultiVector()->Update(alpha, *tmp, beta));
231 }
232
234 virtual bool hasTransposeApply() const {
235 // We do not currently use transpose, try setting it
236 int err = op_->SetUseTranspose(true);
237 op_->SetUseTranspose(false);
238 return (err == 0);
239 }
240
242
244
245
247 std::string description() const { XPETRA_MONITOR("EpetraOperator::description"); return "Epetra_Operator"; }
248
251 XPETRA_MONITOR("EpetraOperator::describe");
252 out << "Epetra_Operator" << std::endl;
253 }
254
256
258
259
261 EpetraInverseOperator(const Teuchos::RCP<Epetra_Operator> &op) : op_(op) { } //TODO removed const
262
263
272
274
275 private:
278
279 };
280
281} // Xpetra namespace
282
283#endif // XPETRA_EPETRAOPERATOR_HPP
#define XPETRA_MONITOR(funcName)
#define XPETRA_ERR_CHECK(arg)
#define XPETRA_DYNAMIC_CAST(type, obj, newObj, exceptionMsg)
static const EVerbosityLevel verbLevel_default
virtual void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const
Computes the operator-multivector application.
virtual Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
The Map associated with the domain of this operator, which must be compatible with X....
std::string description() const
A simple one-line description of this object.
void residual(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &R) const
Compute a residual R = B - (*this) * X.
RCP< Epetra_Operator > op_
The Tpetra::Operator which this class wraps.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with the given verbosity level to a FancyOStream.
virtual Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
The Map associated with the range of this operator, which must be compatible with Y....
EpetraInverseOperator(const Teuchos::RCP< Epetra_Operator > &op)
EpetraOperator constructor to wrap a Epetra_Operator object.
virtual bool hasTransposeApply() const
Whether this operator supports applying the transpose or conjugate transpose.
virtual Teuchos::RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
The Map associated with the domain of this operator, which must be compatible with X....
void residual(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &R) const
Compute a residual R = B - (*this) * X.
std::string description() const
A simple one-line description of this object.
virtual bool hasTransposeApply() const
Whether this operator supports applying the transpose or conjugate transpose.
virtual void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const
Computes the operator-multivector application.
EpetraOperator(const Teuchos::RCP< Epetra_Operator > &op)
EpetraOperator constructor to wrap a Epetra_Operator object.
RCP< Epetra_Operator > op_
The Tpetra::Operator which this class wraps.
EpetraGlobalOrdinal GlobalOrdinal
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with the given verbosity level to a FancyOStream.
virtual Teuchos::RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
The Map associated with the range of this operator, which must be compatible with Y....
Exception throws when you call an unimplemented method of Xpetra.
Exception throws to report errors in the internal logical of the program.
virtual void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)=0
Update multi-vector values with scaled values of A, this = beta*this + alpha*A.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Xpetra namespace