Teko Version of the Day
Loading...
Searching...
No Matches
Teko_SmootherPreconditionerFactory.cpp
1/*
2// @HEADER
3//
4// ***********************************************************************
5//
6// Teko: A package for block and physics based preconditioning
7// Copyright 2010 Sandia Corporation
8//
9// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
10// the U.S. Government retains certain rights in this software.
11//
12// Redistribution and use in source and binary forms, with or without
13// modification, are permitted provided that the following conditions are
14// met:
15//
16// 1. Redistributions of source code must retain the above copyright
17// notice, this list of conditions and the following disclaimer.
18//
19// 2. Redistributions in binary form must reproduce the above copyright
20// notice, this list of conditions and the following disclaimer in the
21// documentation and/or other materials provided with the distribution.
22//
23// 3. Neither the name of the Corporation nor the names of the
24// contributors may be used to endorse or promote products derived from
25// this software without specific prior written permission.
26//
27// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
28// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
29// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
30// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
31// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
32// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
33// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
34// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
35// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
36// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
37// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
38//
39// Questions? Contact Eric C. Cyr (eccyr@sandia.gov)
40//
41// ***********************************************************************
42//
43// @HEADER
44
45*/
46
47// Teko includes
48#include "Teko_SmootherPreconditionerFactory.hpp"
49
50#include "Teko_PreconditionerInverseFactory.hpp"
51#include "Thyra_MultiVectorStdOps.hpp"
52
53namespace Teko {
54
55// A small utility class to help distinguish smoothers
56class SmootherRequestMesg : public RequestMesg {
57public:
58 SmootherRequestMesg(unsigned int block)
59 : RequestMesg("__smoother_request_message__")
60 , block_(block) {}
61
62 unsigned int getBlock() const
63 { return block_; }
64
65private:
66 unsigned int block_;
67};
68
72
73SmootherLinearOp::SmootherLinearOp(const LinearOp & A, const LinearOp & invM,unsigned int applications,bool useDestAsInitialGuess)
74 : A_(A), invM_(invM), applications_(applications), initialGuessType_(Unspecified), requestMesg_(Teuchos::null)
75{
76 // set initial guess type
77 initialGuessType_ = useDestAsInitialGuess ? DestAsInitialGuess : NoInitialGuess;
78}
79
80SmootherLinearOp::SmootherLinearOp(const LinearOp & A, const LinearOp & invM,unsigned int applications,unsigned int block)
81 : A_(A), invM_(invM), applications_(applications), initialGuessType_(RequestInitialGuess), requestMesg_(Teuchos::null)
82{
83 requestMesg_ = Teuchos::rcp(new SmootherRequestMesg(block));
84}
85
86void SmootherLinearOp::implicitApply(const MultiVector & b, MultiVector & x,
87 const double alpha, const double beta) const
88{
89 using Teuchos::RCP;
90
91 MultiVector residual = deepcopy(b); // residual = b
92 MultiVector scrap = deepcopy(b); // scrap = b
93 MultiVector error; // location for initial guess
94
95 // construct initial guess: required to assign starting point for destination
96 // vector appropriately
97 switch(initialGuessType_) {
98 case RequestInitialGuess:
99 // get initial guess from request handler
100 error = deepcopy(getRequestHandler()->request<MultiVector>(*requestMesg_));
101 Thyra::assign<double>(x.ptr(),*error); // x = error (initial guess)
102 break;
103 case DestAsInitialGuess:
104 error = deepcopy(x); // error = x
105 break;
106 case NoInitialGuess:
107 Thyra::assign<double>(x.ptr(),0.0); // x = 0
108 error = deepcopy(x); // error = x
109 break;
110 case Unspecified:
111 default:
112 TEUCHOS_ASSERT(false);
113 }
114
115 for(unsigned int current=0;current<applications_;++current) {
116 // compute current residual
117 Teko::applyOp(A_,error,scrap);
118 Teko::update(-1.0,scrap,1.0,residual); // residual = residual-A*error
119
120 // compute appoximate correction using residual
121 Thyra::assign(error.ptr(),0.0); // set error to zero
122 Teko::applyOp(invM_,residual,error);
123
124 // update solution with error
125 Teko::update(1.0,error,1.0,x); // x = x+error
126 }
127}
128
130void SmootherLinearOp::setRequestHandler(const Teuchos::RCP<RequestHandler> & rh)
131{
132 Teko_DEBUG_SCOPE("SmootherLinearOp::setRequestHandler",10);
133 requestHandler_ = rh;
134}
135
137Teuchos::RCP<RequestHandler> SmootherLinearOp::getRequestHandler() const
138{
139 Teko_DEBUG_SCOPE("SmootherLinearOp::getRequestHandler",10);
140 return requestHandler_;
141}
142
143LinearOp buildSmootherLinearOp(const LinearOp & A,const LinearOp & invM,unsigned int applications,bool useDestAsInitialGuess)
144{
145 return Teuchos::rcp(new SmootherLinearOp(A,invM,applications,useDestAsInitialGuess));
146}
147
148LinearOp buildSmootherLinearOp(const LinearOp & A,const LinearOp & invM,unsigned int applications,unsigned int initialGuessBlock)
149{
150 return Teuchos::rcp(new SmootherLinearOp(A,invM,applications,initialGuessBlock));
151}
152
153
157
159SmootherPreconditionerFactory::SmootherPreconditionerFactory()
160 : sweepCount_(0), initialGuessType_(Unspecified), initialGuessBlock_(0), precFactory_(Teuchos::null)
161{ }
162
166LinearOp SmootherPreconditionerFactory::buildPreconditionerOperator(LinearOp & lo,PreconditionerState & state) const
167{
168 TEUCHOS_TEST_FOR_EXCEPTION(precFactory_==Teuchos::null,std::runtime_error,
169 "ERROR: Teko::SmootherPreconditionerFactory::buildPreconditionerOperator requires that a "
170 << "preconditioner factory has been set. Currently it is null!");
171
172 // preconditions
173 TEUCHOS_ASSERT(sweepCount_>0);
174 TEUCHOS_ASSERT(initialGuessType_!=Unspecified);
175 TEUCHOS_ASSERT(precFactory_!=Teuchos::null);
176
177 // build user specified preconditioner
178 ModifiableLinearOp & invM = state.getModifiableOp("prec");
179 if(invM==Teuchos::null)
180 invM = Teko::buildInverse(*precFactory_,lo);
181 else
182 Teko::rebuildInverse(*precFactory_,lo,invM);
183
184 // conditional on initial guess type, build the smoother
185 switch(initialGuessType_) {
186 case RequestInitialGuess:
187 return buildSmootherLinearOp(lo,invM,sweepCount_,initialGuessBlock_);
188 case DestAsInitialGuess:
189 return buildSmootherLinearOp(lo,invM,sweepCount_,true); // use an initial guess
190 case NoInitialGuess:
191 return buildSmootherLinearOp(lo,invM,sweepCount_,false); // no initial guess
192 case Unspecified:
193 default:
194 TEUCHOS_ASSERT(false);
195 }
196
197 // should never get here
198 TEUCHOS_ASSERT(false);
199 return Teuchos::null;
200}
201
205void SmootherPreconditionerFactory::initializeFromParameterList(const Teuchos::ParameterList & settings)
206{
207 // declare all strings used by this initialization routine
209
210 const std::string str_sweepCount = "Sweep Count";
211 const std::string str_initialGuessBlock = "Initial Guess Block";
212 const std::string str_destAsInitialGuess = "Destination As Initial Guess";
213 const std::string str_precType = "Preconditioner Type";
214
215 // default parameters
217
218 initialGuessType_ = Unspecified;
219 initialGuessBlock_ = 0;
220 sweepCount_ = 0;
221 precFactory_ = Teuchos::null;
222
223 // get sweep counts
225
226 if(settings.isParameter(str_sweepCount))
227 sweepCount_ = settings.get<int>(str_sweepCount);
228
229 // get initial guess
231
232 if(settings.isParameter(str_initialGuessBlock)) {
233 initialGuessBlock_ = settings.get<int>(str_initialGuessBlock);
234 initialGuessType_ = RequestInitialGuess;
235 }
236
237 if(settings.isParameter(str_destAsInitialGuess)) {
238 bool useDest = settings.get<bool>(str_destAsInitialGuess);
239 if(useDest) {
240 TEUCHOS_TEST_FOR_EXCEPTION(initialGuessType_!=Unspecified, std::runtime_error,
241 "Cannot set both \"" << str_initialGuessBlock <<
242 "\" and \"" << str_destAsInitialGuess << "\"");
243
244 initialGuessType_ = DestAsInitialGuess;
245 }
246 }
247
248 // safe to assume if the other values are not turned on there is no initial guess
249 if(initialGuessType_==Unspecified)
250 initialGuessType_ = NoInitialGuess;
251
252 // get preconditioner factory
254
255 TEUCHOS_TEST_FOR_EXCEPTION(not settings.isParameter(str_precType),std::runtime_error,
256 "Parameter \"" << str_precType << "\" is required by a Teko::SmootherPreconditionerFactory");
257
258 // grab library and preconditioner name
259 Teuchos::RCP<const InverseLibrary> il = getInverseLibrary();
260 std::string precName = settings.get<std::string>(str_precType);
261
262 // build preconditioner factory
263 precFactory_ = il->getInverseFactory(precName);
264 TEUCHOS_TEST_FOR_EXCEPTION(precFactory_==Teuchos::null,std::runtime_error,
265 "ERROR: \"" << str_precType << "\" = " << precName
266 << " could not be found");
267
268 // post conditions required to be satisfied
270
271 TEUCHOS_ASSERT(sweepCount_>0);
272 TEUCHOS_ASSERT(initialGuessType_!=Unspecified);
273 TEUCHOS_ASSERT(precFactory_!=Teuchos::null);
274}
275
276} // end namespace Teko
InverseLinearOp buildInverse(const InverseFactory &factory, const LinearOp &A)
Build an inverse operator using a factory and a linear operator.
void rebuildInverse(const InverseFactory &factory, const LinearOp &A, InverseLinearOp &invA)
MultiVector deepcopy(const MultiVector &v)
Perform a deep copy of the vector.
Teuchos::RCP< const InverseLibrary > getInverseLibrary() const
Get the inverse library used by this preconditioner factory.