MueLu Version of the Day
Loading...
Searching...
No Matches
MueLu_MatlabSmoother_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
47#ifndef MUELU_MATLABSMOOTHER_DEF_HPP
48#define MUELU_MATLABSMOOTHER_DEF_HPP
50
51#if defined(HAVE_MUELU_TPETRA) && defined(HAVE_MUELU_MATLAB)
52#include "MueLu_Monitor.hpp"
53
54
55namespace MueLu {
56
57 template <class Scalar,class LocalOrdinal, class GlobalOrdinal, class Node>
59 {
60 SetParameterList(paramList);
61 }
62
63 template <class Scalar,class LocalOrdinal, class GlobalOrdinal, class Node>
65 {
67 ParameterList& pL = const_cast<ParameterList&>(this->GetParameterList());
68 setupFunction_ = pL.get("Setup Function","");
69 solveFunction_ = pL.get("Solve Function","");
70 solveDataSize_ = pL.get("Number of Solver Args", 0);
71 }
72
73 template <class Scalar,class LocalOrdinal, class GlobalOrdinal, class Node>
75 {
76 using namespace std;
77 this->Input(currentLevel, "A");
78 ParameterList& pL = const_cast<ParameterList&>(this->GetParameterList());
79 needsSetup_ = pL.get<string>("Needs");
80 vector<string> needsList = tokenizeList(needsSetup_);
81 for(size_t i = 0; i < needsList.size(); i++)
82 {
83 if(!IsParamMuemexVariable(needsList[i]) && needsList[i] != "Level")
84 this->Input(currentLevel, needsList[i]);
85 }
86 }
87
88 template <class Scalar,class LocalOrdinal, class GlobalOrdinal, class Node>
90 {
91 using namespace std;
92 FactoryMonitor m(*this, "Setup Smoother", currentLevel);
93 if (this->IsSetup() == true)
94 this->GetOStream(Warnings0) << "MueLu::MatlabSmoother::Setup(): Setup() has already been called";
95 vector<RCP<MuemexArg>> InputArgs = processNeeds<Scalar, LocalOrdinal, GlobalOrdinal, Node>(this, needsSetup_, currentLevel);
96 A_ = Factory::Get<RCP<Matrix>>(currentLevel, "A");
97 RCP<MuemexArg> AmatArg = rcp_implicit_cast<MuemexArg>(rcp(new MuemexData<RCP<Matrix>>(A_)));
98 //Always add A to the beginning of InputArgs
99 InputArgs.insert(InputArgs.begin(), AmatArg);
100 // Call mex function
101 if(!setupFunction_.length())
102 throw runtime_error("Invalid matlab function name");
103 solveData_= callMatlab(setupFunction_, solveDataSize_, InputArgs);
104 this->GetOStream(Statistics1) << description() << endl;
105 this->IsSetup(true); //mark the smoother as set up
106 }
107
108 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
109 void MatlabSmoother<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Apply(MultiVector& X, const MultiVector& B, bool InitialGuessIsZero) const
110 {
111 TEUCHOS_TEST_FOR_EXCEPTION(SmootherPrototype::IsSetup() == false, Exceptions::RuntimeError,
112 "MueLu::MatlabSmoother::Apply(): Setup() has not been called");
113 using namespace Teuchos;
114 using namespace std;
115 if(InitialGuessIsZero)
116 X.putScalar(0.0);
117 // Push on A as first input
118 vector<RCP<MuemexArg>> InputArgs;
119 InputArgs.push_back(rcp(new MuemexData<RCP<Matrix>>(A_)));
120 // Push on LHS & RHS
121 RCP<MultiVector> Xrcp(&X, false);
122 MultiVector* BPtrNonConst = (MultiVector*) &B;
123 RCP<MultiVector> Brcp = rcp<MultiVector>(BPtrNonConst, false);
124 RCP<MuemexData<RCP<MultiVector>>> XData = rcp(new MuemexData<RCP<MultiVector>>(Xrcp));
125 RCP<MuemexData<RCP<MultiVector>>> BData = rcp(new MuemexData<RCP<MultiVector>>(Brcp));
126 InputArgs.push_back(XData);
127 InputArgs.push_back(BData);
128 for(size_t i = 0; i < solveData_.size(); i++)
129 InputArgs.push_back(solveData_[i]);
130 if(!solveFunction_.length()) throw std::runtime_error("Invalid matlab function name");
131 vector<Teuchos::RCP<MuemexArg>> mexOutput = callMatlab(solveFunction_, 1, InputArgs);
132 RCP<MuemexData<RCP<MultiVector>>> mydata = Teuchos::rcp_static_cast<MuemexData<RCP<MultiVector>>>(mexOutput[0]);
133 X = *(mydata->getData());
134 }
135
136 template <class Scalar,class LocalOrdinal, class GlobalOrdinal, class Node>
137 RCP<MueLu::SmootherPrototype<Scalar, LocalOrdinal, GlobalOrdinal, Node>> MatlabSmoother<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Copy() const
138 {
139 RCP<MatlabSmoother> smoother = rcp(new MatlabSmoother(*this) );
140 smoother->SetParameterList(this->GetParameterList());
141 return smoother;
142 }
143
144 template <class Scalar,class LocalOrdinal, class GlobalOrdinal, class Node>
146 std::ostringstream out;
148 out << "Matlab Smoother("<<setupFunction_<<"/"<<solveFunction_<<")";
149 } else {
151 }
152 return out.str();
153 }
154
155 template <class Scalar,class LocalOrdinal, class GlobalOrdinal, class Node>
156 void MatlabSmoother<Scalar, LocalOrdinal, GlobalOrdinal, Node>::print(Teuchos::FancyOStream &out, const VerbLevel verbLevel) const {
158
159 if (verbLevel & Parameters0)
160 out << "Matlab Smoother("<<setupFunction_<<"/"<<solveFunction_<<")";
161
162 if (verbLevel & Parameters1) {
163 out0 << "Parameter list: " << std::endl;
164 Teuchos::OSTab tab2(out);
165 out << this->GetParameterList();
166 }
167
168 if (verbLevel & Debug) {
169 out0 << "IsSetup: " << Teuchos::toString(SmootherPrototype::IsSetup()) << std::endl;
170 }
171 }
172
173
174// Dummy specializations for GO = long long
175/*template <>
176void MatlabSmoother<double,int,long long>::Setup(Level& currentLevel) {
177 throw std::runtime_error("MatlabSmoother does not support GlobalOrdinal == long long.");
178}
179template <>
180void MatlabSmoother<std::complex<double>,int,long long>::Setup(Level& currentLevel) {
181 throw std::runtime_error("MatlabSmoother does not support GlobalOrdinal == long long.");
182}
183
184template <>
185void MatlabSmoother<double,int,long long>::Apply(MultiVector& X, const MultiVector& B, bool InitialGuessIsZero) const {
186 throw std::runtime_error("MatlabSmoother does not support GlobalOrdinal == long long.");
187}
188template <>
189void MatlabSmoother<std::complex<double>,int,long long>::Apply(MultiVector& X, const MultiVector& B, bool InitialGuessIsZero) const {
190 throw std::runtime_error("MatlabSmoother does not support GlobalOrdinal == long long.");
191}*/
192
193
194} // namespace MueLu
195
196#endif // HAVE_MUELU_TPETRA && HAVE_MUELU_MATLAB
197#endif // MUELU_MATLABSMOOTHER_DEF_HPP
#define MUELU_DESCRIBE
Helper macro for implementing Describable::describe() for BaseClass objects.
virtual std::string description() const
Return a simple one-line description of this object.
Exception throws to report errors in the internal logical of the program.
Timer to be used in factories. Similar to Monitor but with additional timers.
Class that holds all level-specific information.
Definition: MueLu_Level.hpp:99
void DeclareInput(Level &currentLevel) const
Input.
void print(Teuchos::FancyOStream &out, const VerbLevel verbLevel=Default) const
Print the object with some verbosity level to an FancyOStream object.
void Setup(Level &currentLevel)
Set up the smoother.
RCP< SmootherPrototype > Copy() const
friend class MatlabSmoother
Constructor.
void Apply(MultiVector &X, const MultiVector &B, bool InitialGuessIsZero=false) const
Apply the preconditioner.
void SetParameterList(const Teuchos::ParameterList &paramList)
Set parameters from a parameter list and return with default values.
std::string description() const
Return a simple one-line description of this object.
virtual void SetParameterList(const Teuchos::ParameterList &paramList)=0
Set parameters from a parameter list and return with default values.
bool IsSetup() const
Get the state of a smoother prototype.
Namespace for MueLu classes and methods.
bool IsParamMuemexVariable(const std::string &name)
@ Warnings0
Important warning messages (one line)
@ Debug
Print additional debugging information.
@ Statistics1
Print more statistics.
@ Parameters0
Print class parameters.
@ Parameters1
Print class parameters (more parameters, more verbose)
std::vector< RCP< MuemexArg > > callMatlab(std::string function, int numOutputs, std::vector< RCP< MuemexArg > > args)
std::vector< std::string > tokenizeList(const std::string &params)