MueLu Version of the Day
Loading...
Searching...
No Matches
MueLu_TogglePFactory_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#ifndef MUELU_TOGGLEPFACTORY_DEF_HPP
47#define MUELU_TOGGLEPFACTORY_DEF_HPP
48
49#include <Xpetra_Matrix.hpp>
50
52
54#include "MueLu_Level.hpp"
55#include "MueLu_MasterList.hpp"
56#include "MueLu_Monitor.hpp"
58
59namespace MueLu {
60
61 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
63 RCP<ParameterList> validParamList = rcp(new ParameterList());
64
65#define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name))
66 SET_VALID_ENTRY("toggle: mode");
67 SET_VALID_ENTRY("semicoarsen: number of levels");
68#undef SET_VALID_ENTRY
69
70 return validParamList;
71 }
72
73 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
75 // request/release "P" and coarse level "Nullspace"
76 for (std::vector<RCP<const FactoryBase> >::const_iterator it = prolongatorFacts_.begin(); it != prolongatorFacts_.end(); ++it) {
77 coarseLevel.DeclareInput("P", (*it).get(), this); // request/release "P" (dependencies are not affected)
78 coarseLevel.DeclareInput("RfromPfactory", (*it).get(), this);
79 (*it)->CallDeclareInput(coarseLevel); // request dependencies
80 }
81 for (std::vector<RCP<const FactoryBase> >::const_iterator it = ptentFacts_.begin(); it != ptentFacts_.end(); ++it) {
82 coarseLevel.DeclareInput("P", (*it).get(), this); // request/release "Ptent" (dependencies are not affected)
83 (*it)->CallDeclareInput(coarseLevel); // request dependencies
84 }
85 for (std::vector<RCP<const FactoryBase> >::const_iterator it = nspFacts_.begin(); it != nspFacts_.end(); ++it) {
86 coarseLevel.DeclareInput("Nullspace", (*it).get(), this); // request/release coarse "Nullspace" (dependencies are not affected)
87 (*it)->CallDeclareInput(coarseLevel); // request dependencies
88 }
89
90 // The factory needs the information about the number of z-layers. While this information is
91 // provided by the user for the finest level, the factory itself is responsible to provide the
92 // corresponding information on the coarser levels. Since a factory cannot be dependent on itself
93 // we use the NoFactory class as generator class, but remove the UserData keep flag, such that
94 // "NumZLayers" is part of the request/release mechanism.
95 // Please note, that this prevents us from having several (independent) CoarsePFactory instances!
96 // TODO: allow factory to dependent on self-generated data for TwoLevelFactories -> introduce ExpertRequest/Release in Level
97 fineLevel.DeclareInput("NumZLayers", NoFactory::get(), this);
98 fineLevel.RemoveKeepFlag("NumZLayers", NoFactory::get(), MueLu::UserData);
99
100 hasDeclaredInput_ = true;
101 }
102
103 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
105 FactoryMonitor m(*this, "Prolongator toggle", coarseLevel);
106 std::ostringstream levelstr;
107 levelstr << coarseLevel.GetLevelID();
108
109 TEUCHOS_TEST_FOR_EXCEPTION(nspFacts_.size() != prolongatorFacts_.size(), Exceptions::RuntimeError, "MueLu::TogglePFactory::Build: The number of provided prolongator factories and coarse nullspace factories must be identical.");
110 TEUCHOS_TEST_FOR_EXCEPTION(nspFacts_.size() != 2, Exceptions::RuntimeError, "MueLu::TogglePFactory::Build: TogglePFactory needs two different transfer operator strategies for toggling."); // TODO adapt this/weaken this as soon as other toggling strategies are introduced.
111
112 // decision routine which prolongator factory to be used
113 int nProlongatorFactory = 0; // default behavior: use first prolongator in list
114
115 // extract user parameters
116 const Teuchos::ParameterList & pL = GetParameterList();
117 std::string mode = Teuchos::as<std::string>(pL.get<std::string>("toggle: mode"));
118 int semicoarsen_levels = Teuchos::as<int>(pL.get<int>("semicoarsen: number of levels"));
119
120 TEUCHOS_TEST_FOR_EXCEPTION(mode!="semicoarsen", Exceptions::RuntimeError, "MueLu::TogglePFactory::Build: The 'toggle: mode' parameter must be set to 'semicoarsen'. No other mode supported, yet.");
121
122 LO NumZDir = -1;
123 if(fineLevel.IsAvailable("NumZLayers", NoFactory::get())) {
124 NumZDir = fineLevel.Get<LO>("NumZLayers", NoFactory::get()); //obtain info
125 GetOStream(Runtime1) << "Number of layers for semicoarsening: " << NumZDir << std::endl;
126 }
127
128 // Make a decision which prolongator to be used.
129 if(fineLevel.GetLevelID() >= semicoarsen_levels || NumZDir == 1) {
130 nProlongatorFactory = 1;
131 } else {
132 nProlongatorFactory = 0;
133 }
134
135 RCP<Matrix> P = Teuchos::null;
136 RCP<Matrix> Ptent = Teuchos::null;
137 RCP<MultiVector> coarseNullspace = Teuchos::null;
138
139 // call Build for selected transfer operator
140 GetOStream(Runtime0) << "TogglePFactory: call transfer factory: " << (prolongatorFacts_[nProlongatorFactory])->description() << std::endl;
141 prolongatorFacts_[nProlongatorFactory]->CallBuild(coarseLevel);
142 P = coarseLevel.Get< RCP<Matrix> >("P", (prolongatorFacts_[nProlongatorFactory]).get());
143 RCP<Matrix> R = Teuchos::null;
144 int Rplaceholder = -1; // Used to indicate that an R matrix has not been produced by a prolongator factory, but
145 // that it is capable of producing one and should be invoked a 2nd time in restrictor mode
146 // (e.g. with PgPFactory). prolongatorFacts_[Rplaceholder] is factory that can produce R
147 // matrix, which might be later invoked by MueLu_RfromP_Or_TransP
148 if (coarseLevel.IsAvailable("RfromPfactory",(prolongatorFacts_[nProlongatorFactory]).get())) {
149 std::string strType = coarseLevel.GetTypeName("RfromPfactory", (prolongatorFacts_[nProlongatorFactory]).get());
150 if (strType == "int") Rplaceholder = nProlongatorFactory;
151 else R = coarseLevel.Get< RCP<Matrix> >("RfromPfactory", (prolongatorFacts_[nProlongatorFactory]).get());
152 // Need to get R (and set it below) so that TogglePFactory is given credit for creating R
153 }
154 // do not call "Build" for "Ptent" factory since it should automatically be called recursively
155 // through the "Build" call for "P"
156 Ptent = coarseLevel.Get< RCP<Matrix> >("P", (ptentFacts_[nProlongatorFactory]).get());
157 coarseNullspace = coarseLevel.Get< RCP<MultiVector> >("Nullspace", (nspFacts_[nProlongatorFactory]).get());
158
159 // Release dependencies of all prolongator and coarse level null spaces
160 for(size_t t=0; t<nspFacts_.size(); ++t) {
161 coarseLevel.Release(*(prolongatorFacts_[t]));
162 coarseLevel.Release(*(ptentFacts_[t]));
163 coarseLevel.Release(*(nspFacts_[t]));
164 }
165
166 // store prolongator with this factory identification.
167 Set(coarseLevel, "P", P);
168 // Three cases:
169 // 1) R already computed and TogglePFactory takes credit for constructing it
170 // 2) R not computed but prolongatorFacts_[Rplaceholder] can produce it
171 // 3) R not computed and prolongator can not produce it
172 if (R != Teuchos::null) Set(coarseLevel, "RfromPfactory", R);
173 else if (Rplaceholder != -1) Set(coarseLevel, "RfromPfactory", Teuchos::as<int>(Rplaceholder));
174 Set(coarseLevel, "Nullspace", coarseNullspace);
175 Set(coarseLevel, "Ptent", Ptent);
176 Set(coarseLevel, "Chosen P", nProlongatorFactory);
177 } //Build()
178
179 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
181 // check if it's a TwoLevelFactoryBase based transfer factory
182 TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::rcp_dynamic_cast<const TwoLevelFactoryBase>(factory) == Teuchos::null, Exceptions::BadCast,
183 "MueLu::TogglePFactory::AddProlongatorFactory: Transfer factory is not derived from TwoLevelFactoryBase. "
184 "This is very strange. (Note: you can remove this exception if there's a good reason for)");
185 TEUCHOS_TEST_FOR_EXCEPTION(hasDeclaredInput_, Exceptions::RuntimeError, "MueLu::TogglePFactory::AddProlongatorFactory: Factory is being added after we have already declared input");
186 prolongatorFacts_.push_back(factory);
187 }
188
189 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
191 // check if it's a TwoLevelFactoryBase based transfer factory
192 TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::rcp_dynamic_cast<const TwoLevelFactoryBase>(factory) == Teuchos::null, Exceptions::BadCast,
193 "MueLu::TogglePFactory::AddPtentFactory: Transfer factory is not derived from TwoLevelFactoryBase. "
194 "This is very strange. (Note: you can remove this exception if there's a good reason for)");
195 TEUCHOS_TEST_FOR_EXCEPTION(hasDeclaredInput_, Exceptions::RuntimeError, "MueLu::TogglePFactory::AddPtentFactory: Factory is being added after we have already declared input");
196 ptentFacts_.push_back(factory);
197 }
198
199 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
201 // check if it's a TwoLevelFactoryBase based transfer factory
202 TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::rcp_dynamic_cast<const TwoLevelFactoryBase>(factory) == Teuchos::null, Exceptions::BadCast,
203 "MueLu::TogglePFactory::AddCoarseNullspaceFactory: Transfer factory is not derived from TwoLevelFactoryBase. Make sure you provide the factory which generates the coarse level nullspace information. Usually this is a prolongator factory."
204 "This is very strange. (Note: you can remove this exception if there's a good reason for)");
205 TEUCHOS_TEST_FOR_EXCEPTION(hasDeclaredInput_, Exceptions::RuntimeError, "MueLu::TogglePFactory::AddCoarseNullspaceFactory: Factory is being added after we have already declared input");
206 nspFacts_.push_back(factory);
207 }
208
209
210} //namespace MueLu
211
212#endif // MUELU_TOGGLEPFACTORY_DEF_HPP
#define SET_VALID_ENTRY(name)
Exception indicating invalid cast attempted.
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.
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()
void Release(const FactoryBase &factory)
Decrement the storage counter for all the inputs of a factory.
void RemoveKeepFlag(const std::string &ename, const FactoryBase *factory, KeepType keep=MueLu::All)
int GetLevelID() const
Return level number.
std::string GetTypeName(const std::string &ename, const FactoryBase *factory=NoFactory::get())
GetTypeName returns type string of variable stored using ename and factory.
T & Get(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Get data without decrementing associated storage counter (i.e., read-only access)....
static const NoFactory * get()
void AddCoarseNullspaceFactory(const RCP< const FactoryBase > &factory)
Add a coarse nullspace factory in the end of list of coarse nullspace factories.
void AddProlongatorFactory(const RCP< const FactoryBase > &factory)
Add a prolongator factory in the end of list of prolongator factories.
void AddPtentFactory(const RCP< const FactoryBase > &factory)
Add a tentative prolongator factory in the end of list of prolongator factories.
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
void Build(Level &fineLevel, Level &coarseLevel) const
Build method.
Namespace for MueLu classes and methods.
@ UserData
User data are always kept. This flag is set automatically when Level::Set("data", data) is used....
@ Runtime0
One-liner description of what is happening.
@ Runtime1
Description of what is happening (more verbose)