MueLu Version of the Day
Loading...
Searching...
No Matches
MueLu_AMGXOperator_def.hpp
Go to the documentation of this file.
1// @HEADER
2//
3// ***********************************************************************
4//
5// MueLu: A package for multigrid based preconditioning
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
47#ifndef MUELU_AMGXOPERATOR_DEF_HPP
48#define MUELU_AMGXOPERATOR_DEF_HPP
49
50
51#if defined (HAVE_MUELU_AMGX)
53
54namespace MueLu {
55
56 template<class Node>
57 Teuchos::RCP<const Tpetra::Map<int,int,Node> >
59 return domainMap_;
60 }
61
62 template<class Node>
63 Teuchos::RCP<const Tpetra::Map<int,int,Node> > AMGXOperator<double,int,int,Node>::getRangeMap() const {
64 return rangeMap_;
65 }
66
67 template<class Node>
68 void AMGXOperator<double,int,int,Node>::apply(const Tpetra::MultiVector<double,int,int,Node>& X,
69 Tpetra::MultiVector<double,int,int,Node>& Y,
70 Teuchos::ETransp mode, double alpha, double beta) const {
71
72 RCP<const Teuchos::Comm<int> > comm = Y.getMap()->getComm();
73
74 ArrayRCP<const double> mueluXdata, amgxXdata;
75 ArrayRCP<double> mueluYdata, amgxYdata;
76
77 try {
78 for (int i = 0; i < (int)Y.getNumVectors(); i++) {
79 {
80 vectorTimer1_->start();
81
82 mueluXdata = X.getData(i);
83 mueluYdata = Y.getDataNonConst(i);
84
85 if (comm->getSize() == 1) {
86 amgxXdata = mueluXdata;
87 amgxYdata = mueluYdata;
88
89 } else {
90 int n = mueluXdata.size();
91
92 amgxXdata.resize(n);
93 amgxYdata.resize(n);
94
95 ArrayRCP<double> amgxXdata_nonConst = Teuchos::arcp_const_cast<double>(amgxXdata);
96 for (int j = 0; j < n; j++) {
97 amgxXdata_nonConst[muelu2amgx_[j]] = mueluXdata[j];
98 amgxYdata [muelu2amgx_[j]] = mueluYdata[j];
99 }
100 }
101
102 AMGX_vector_upload(X_, N_, 1, &amgxXdata[0]);
103 AMGX_vector_upload(Y_, N_, 1, &amgxYdata[0]);
104
105 vectorTimer1_->stop();
106 vectorTimer1_->incrementNumCalls();
107 }
108
109 // Solve the system and time.
110 solverTimer_->start();
111 AMGX_solver_solve(Solver_, X_, Y_);
112 solverTimer_->stop();
113 solverTimer_->incrementNumCalls();
114
115 {
116 vectorTimer2_->start();
117
118 AMGX_vector_download(Y_, &amgxYdata[0]);
119
120 if (comm->getSize() > 1) {
121 int n = mueluYdata.size();
122
123 for (int j = 0; j < n; j++)
124 mueluYdata[j] = amgxYdata[muelu2amgx_[j]];
125 }
126
127 vectorTimer2_->stop();
128 vectorTimer2_->incrementNumCalls();
129 }
130 }
131
132 } catch (std::exception& e) {
133 std::string errMsg = std::string("Caught an exception in MueLu::AMGXOperator::Apply():\n") + e.what() + "\n";
134 throw Exceptions::RuntimeError(errMsg);
135 }
136 }
137
138 template<class Node>
140 return false;
141 }
142
143} // namespace
144#endif //if defined(HAVE_MUELU_AMGX)
145
146#endif //ifdef MUELU_AMGXOPERATOR_DEF_HPP
void apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const
Returns a solution for the linear system AX=Y in the Tpetra::MultiVector X.
Teuchos::RCP< const Map > getRangeMap() const
Returns the Tpetra::Map object associated with the range of this operator.
bool hasTransposeApply() const
Indicates whether this operator supports applying the adjoint operator.
Teuchos::RCP< const Map > getDomainMap() const
Returns the Tpetra::Map object associated with the domain of this operator.
Exception throws to report errors in the internal logical of the program.
Namespace for MueLu classes and methods.