MueLu  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MueLu_RebalanceBlockAcFactory_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_REBALANCEBLOCKACFACTORY_DEF_HPP_
48 #define MUELU_REBALANCEBLOCKACFACTORY_DEF_HPP_
49 
50 #include <Xpetra_Matrix.hpp>
51 #include <Xpetra_CrsMatrix.hpp>
52 #include <Xpetra_CrsMatrixWrap.hpp>
53 #include <Xpetra_MatrixFactory.hpp>
54 #include <Xpetra_MapExtractor.hpp>
56 #include <Xpetra_StridedMap.hpp>
59 
60 #include <Xpetra_VectorFactory.hpp>
61 
63 
65 #include "MueLu_HierarchyUtils.hpp"
66 #include "MueLu_MasterList.hpp"
67 #include "MueLu_Monitor.hpp"
68 #include "MueLu_PerfUtils.hpp"
69 #include "MueLu_RAPFactory.hpp"
70 
71 namespace MueLu {
72 
73  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
75 
76  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
78  RCP<ParameterList> validParamList = rcp(new ParameterList());
79 
80 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name))
81  SET_VALID_ENTRY("repartition: use subcommunicators");
82 #undef SET_VALID_ENTRY
83 
84  validParamList->set<RCP<const FactoryBase> >("A", Teuchos::null, "Generating factory of the matrix A for rebalancing");
85  validParamList->set<RCP<const FactoryBase> >("Importer", Teuchos::null, "Generating factory of the matrix Importer for rebalancing");
86  validParamList->set<RCP<const FactoryBase> >("SubImporters", Teuchos::null, "Generating factory of the matrix sub-Importers for rebalancing");
87 
88  return validParamList;
89  }
90 
91  template <class Scalar,class LocalOrdinal, class GlobalOrdinal, class Node>
93  FactManager_.push_back(FactManager);
94  }
95 
96  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
98  Input(coarseLevel, "A");
99 
100  std::vector<Teuchos::RCP<const FactoryManagerBase> >::const_iterator it;
101  for(it = FactManager_.begin(); it!=FactManager_.end(); ++it) {
102  SetFactoryManager fineSFM (rcpFromRef(fineLevel), *it);
103  SetFactoryManager coarseSFM(rcpFromRef(coarseLevel), *it);
104 
105  coarseLevel.DeclareInput("Importer",(*it)->GetFactory("Importer").get(), this);
106  }
107 
108  // Use the non-manager path if the maps / importers are generated in one place
109  if(FactManager_.size() == 0) {
110  Input(coarseLevel,"SubImporters");
111  }
112 
113  }
114 
115  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
117  FactoryMonitor m(*this, "Computing blocked Ac", coarseLevel);
118 
119  RCP<Teuchos::FancyOStream> out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
120 
121  RCP<Matrix> originalAc = Get< RCP<Matrix> >(coarseLevel, "A");
122 
124  TEUCHOS_TEST_FOR_EXCEPTION(bA==Teuchos::null, Exceptions::BadCast, "MueLu::RebalanceBlockAcFactory::Build: input matrix A is not of type BlockedCrsMatrix! error.");
125  TEUCHOS_TEST_FOR_EXCEPTION(bA->Rows() != bA->Cols(), Exceptions::RuntimeError, "MueLu::RebalanceBlockAcFactory::Build: Blocked operator has " << bA->Rows() << " and " << bA->Cols() << ". We only support square matrices (with same number of blocks and columns).");
126 
127  // Variables to set up map extractors for blocked operators
128  std::vector<GO> fullRangeMapVector;
129  std::vector<GO> fullDomainMapVector;
130  std::vector<RCP<const Map> > subBlockARangeMaps;
131  std::vector<RCP<const Map> > subBlockADomainMaps;
132  subBlockARangeMaps.reserve(bA->Rows());
133  subBlockADomainMaps.reserve(bA->Cols());
134 
135  // store map extractors
136  Teuchos::RCP<const MapExtractorClass> rangeMapExtractor = bA->getRangeMapExtractor();
137  Teuchos::RCP<const MapExtractorClass> domainMapExtractor = bA->getDomainMapExtractor();
138 
139  // check if GIDs for full maps have to be sorted:
140  // For the Thyra mode ordering they do not have to be sorted since the GIDs are
141  // numbered as 0...n1,0...,n2 (starting with zero for each subblock). The MapExtractor
142  // generates unique GIDs during the construction.
143  // For Xpetra style, the GIDs have to be reordered. Such that one obtains a ordered
144  // list of GIDs in an increasing ordering. In Xpetra, the GIDs are all unique through
145  // out all submaps.
146  bool bThyraRangeGIDs = rangeMapExtractor->getThyraMode();
147  bool bThyraDomainGIDs = domainMapExtractor->getThyraMode();
148 
149  // vector containing rebalanced blocks (for final output)
150  std::vector<RCP<Matrix> > subBlockRebA =
151  std::vector<RCP<Matrix> >(bA->Cols() * bA->Rows(), Teuchos::null);
152 
153  // vector with Import objects from the different
154  // RepartitionFactory instances
155  std::vector<RCP<const Import> > importers = std::vector<RCP<const Import> >(bA->Rows(), Teuchos::null);
156  std::vector<Teuchos::RCP<const FactoryManagerBase> >::const_iterator it;
157  size_t idx = 0;
158  for(it = FactManager_.begin(); it!=FactManager_.end(); ++it) {
159  SetFactoryManager fineSFM (rcpFromRef(fineLevel), *it);
160  SetFactoryManager coarseSFM(rcpFromRef(coarseLevel), *it);
161 
162  RCP<const Import> rebalanceImporter = coarseLevel.Get<RCP<const Import> >("Importer", (*it)->GetFactory("Importer").get());
163  importers[idx] = rebalanceImporter;
164  idx++;
165  }
166  if(FactManager_.size() == 0) {
167  importers = Get<std::vector<RCP<const Import> > >(coarseLevel,"SubImporters");
168  }
169 
170 
171  // restrict communicator?
172  bool bRestrictComm = false;
173  const ParameterList& pL = GetParameterList();
174  if (pL.get<bool>("repartition: use subcommunicators") == true)
175  bRestrictComm = true;
176 
177  RCP<ParameterList> XpetraList = Teuchos::rcp(new ParameterList());
178  if(bRestrictComm)
179  XpetraList->set("Restrict Communicator",true);
180  else
181  XpetraList->set("Restrict Communicator",false);
182 
183  // communicator for final (rebalanced) operator.
184  // If communicator is not restricted it should be identical to the communicator in bA
185  RCP<const Teuchos::Comm<int> > rebalancedComm = Teuchos::null;
186 
187  // loop through all blocks and rebalance blocks
188  // Note: so far we do not support rebalancing of nested operators
189  // TODO add a check for this
190  for(size_t i=0; i<bA->Rows(); i++) {
191  for(size_t j=0; j<bA->Cols(); j++) {
192  // extract matrix block
193  RCP<Matrix> Aij = bA->getMatrix(i, j);
194 
195  std::stringstream ss; ss << "Rebalancing matrix block A(" << i << "," << j << ")";
196  SubFactoryMonitor subM(*this, ss.str(), coarseLevel);
197 
198  RCP<Matrix> rebAij = Teuchos::null;
199  // General rebalancing
200  if( importers[i] != Teuchos::null &&
201  importers[j] != Teuchos::null &&
202  Aij != Teuchos::null) {
203  RCP<const Map> targetRangeMap = importers[i]->getTargetMap();
204  RCP<const Map> targetDomainMap = importers[j]->getTargetMap();
205 
206  // Copy the block Aij
207  // TAW: Do we need a copy or can we do in-place rebalancing?
208  // If we do in-place rebalancing the original distribution is lost
209  // We don't really need it any more, though.
210  //RCP<Matrix> cAij = MatrixFactory::BuildCopy(Aij);
211  RCP<Matrix> cAij = Aij; // do not copy the matrix data (just an rcp pointer)
212 
213  // create a new importer for column map needed for rebalanced Aij
214  Teuchos::RCP<const Import> rebAijImport = ImportFactory::Build(importers[j]->getTargetMap(),cAij->getColMap());
215  TEUCHOS_TEST_FOR_EXCEPTION(rebAijImport.is_null() == true,Exceptions::RuntimeError, "MueLu::RebalanceBlockAcFactory::Build: Importer associated with block " << j << " is null.");
216 
217  Teuchos::RCP<const CrsMatrixWrap> cAwij = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(cAij);
218  TEUCHOS_TEST_FOR_EXCEPTION(cAwij.is_null() == true,Exceptions::RuntimeError, "MueLu::RebalanceBlockAcFactory::Build: Block (" << i << "," << j << ") is not of type CrsMatrix. We cannot rebalanced (nested) operators.");
219  Teuchos::RCP<CrsMatrix> cAmij = cAwij->getCrsMatrix();
220 
221  // change domain map to rebalanced domain map (in-place). Update the importer to represent the column map
222  //cAmij->replaceDomainMapAndImporter(importers[j]->getTargetMap(),rebAijImport);
223 
224  // rebalance rows of matrix block. Don't change the domain map (-> Teuchos::null)
225  // NOTE: If the communicator is restricted away, Build returns Teuchos::null.
226  rebAij = MatrixFactory::Build(cAij, *(importers[i]), *(importers[j]), targetDomainMap, targetRangeMap, XpetraList);
227  } // rebalance matrix block A(i,i)
228  else {
229  rebAij = Aij; // no rebalancing or empty block!
230  }
231 
232  // store new block in output
233  subBlockRebA[i*bA->Cols() + j] = rebAij;
234 
235  if (!rebAij.is_null()) {
236  // store communicator
237  if(rebalancedComm.is_null()) rebalancedComm = rebAij->getRowMap()->getComm();
238 
239  // printout rebalancing information
240  RCP<ParameterList> params = rcp(new ParameterList());
241  params->set("printLoadBalancingInfo", true);
242  std::stringstream ss2; ss2 << "A(" << i << "," << j << ") rebalanced:";
243  GetOStream(Statistics0) << PerfUtils::PrintMatrixInfo(*rebAij, ss2.str(), params);
244  }
245  } // loop over columns j
246 
247  // fix striding information of diagonal blocks
248  // Note: we do not care about the off-diagonal blocks. We just make sure, that the
249  // diagonal blocks have the corresponding striding information from the map extractors
250  // Note: the diagonal block never should be zero.
251  // TODO what if a diagonal block is Teuchos::null?
252  if ( subBlockRebA[i*bA->Cols() + i].is_null() == false ) {
253  RCP<Matrix> rebAii = subBlockRebA[i*bA->Cols() + i];
254  Teuchos::RCP<const StridedMap> orig_stridedRgMap = Teuchos::rcp_dynamic_cast<const StridedMap>(rangeMapExtractor->getMap(i,rangeMapExtractor->getThyraMode()));
255  Teuchos::RCP<const Map> stridedRgMap = Teuchos::null;
256  if(orig_stridedRgMap != Teuchos::null) {
257  std::vector<size_t> stridingData = orig_stridedRgMap->getStridingData();
258  Teuchos::ArrayView< const GlobalOrdinal > nodeRangeMapii = rebAii->getRangeMap()->getNodeElementList();
259  stridedRgMap = StridedMapFactory::Build(
260  bA->getRangeMap()->lib(),
262  nodeRangeMapii,
263  rebAii->getRangeMap()->getIndexBase(),
264  stridingData,
265  rebalancedComm, /*rebAii->getRangeMap()->getComm(),*/ /* restricted communicator */
266  orig_stridedRgMap->getStridedBlockId(),
267  orig_stridedRgMap->getOffset());
268  }
269  Teuchos::RCP<const StridedMap> orig_stridedDoMap = Teuchos::rcp_dynamic_cast<const StridedMap>(domainMapExtractor->getMap(i,domainMapExtractor->getThyraMode()));
270  Teuchos::RCP<const Map> stridedDoMap = Teuchos::null;
271  if(orig_stridedDoMap != Teuchos::null) {
272  std::vector<size_t> stridingData = orig_stridedDoMap->getStridingData();
273  Teuchos::ArrayView< const GlobalOrdinal > nodeDomainMapii = rebAii->getDomainMap()->getNodeElementList();
274  stridedDoMap = StridedMapFactory::Build(
275  bA->getDomainMap()->lib(),
277  nodeDomainMapii,
278  rebAii->getDomainMap()->getIndexBase(),
279  stridingData,
280  rebalancedComm, /*rebAii->getDomainMap()->getComm(), *//* restricted communicator */
281  orig_stridedDoMap->getStridedBlockId(),
282  orig_stridedDoMap->getOffset());
283  }
284 
285  if(bRestrictComm) {
286  stridedRgMap->removeEmptyProcesses();
287  stridedDoMap->removeEmptyProcesses();
288  }
289 
290  TEUCHOS_TEST_FOR_EXCEPTION(stridedRgMap == Teuchos::null,Exceptions::RuntimeError, "MueLu::RebalanceBlockAcFactory::Build: failed to generate striding information. error.");
291  TEUCHOS_TEST_FOR_EXCEPTION(stridedDoMap == Teuchos::null,Exceptions::RuntimeError, "MueLu::RebalanceBlockAcFactory::Build: failed to generate striding information. error.");
292 
293  // replace stridedMaps view in diagonal sub block
294  if(rebAii->IsView("stridedMaps")) rebAii->RemoveView("stridedMaps");
295  rebAii->CreateView("stridedMaps", stridedRgMap, stridedDoMap);
296  // collect Xpetra-based global row ids for map extractors
297  subBlockARangeMaps.push_back(rebAii->getRowMap("stridedMaps"));
298  Teuchos::ArrayView< const GlobalOrdinal > nodeRangeMap = rebAii->getRangeMap()->getNodeElementList();
299  // append the GIDs in the end. Do not sort if we have Thyra style GIDs
300  fullRangeMapVector.insert(fullRangeMapVector.end(), nodeRangeMap.begin(), nodeRangeMap.end());
301  if(bThyraRangeGIDs == false)
302  sort(fullRangeMapVector.begin(), fullRangeMapVector.end());
303 
304  subBlockADomainMaps.push_back(rebAii->getColMap("stridedMaps"));
305  Teuchos::ArrayView< const GlobalOrdinal > nodeDomainMap = rebAii->getDomainMap()->getNodeElementList();
306  // append the GIDs in the end. Do not sort if we have Thyra style GIDs
307  fullDomainMapVector.insert(fullDomainMapVector.end(), nodeDomainMap.begin(), nodeDomainMap.end());
308  if(bThyraDomainGIDs == false)
309  sort(fullDomainMapVector.begin(), fullDomainMapVector.end());
310  } // end if rebAii != Teuchos::null
311  } // loop over rows i
312 
313  // all sub blocks are rebalanced (if available)
314 
315  // Short cut if this processor is not in the list of active processors
316  if (rebalancedComm == Teuchos::null) {
317  GetOStream(Debug,-1) << "RebalanceBlockedAc: deactivate proc " << originalAc->getRowMap()->getComm()->getRank() << std::endl;
318  // TAW: it is important that we create a dummy object of type BlockedCrsMatrix (even if we set it to Teuchos::null)
319  Teuchos::RCP<BlockedCrsMatrix> reb_bA = Teuchos::null;
320  coarseLevel.Set("A", Teuchos::rcp_dynamic_cast<Matrix>(reb_bA), this);
321  return;
322  }
323 
324  // now, subBlockRebA contains all rebalanced matrix blocks
325  // extract map index base from maps of blocked A
326  GO rangeIndexBase = bA->getRangeMap()->getIndexBase();
327  GO domainIndexBase = bA->getDomainMap()->getIndexBase();
328 
329  Teuchos::ArrayView<GO> fullRangeMapGIDs(fullRangeMapVector.size() ? &fullRangeMapVector[0] : 0,fullRangeMapVector.size());
330  Teuchos::RCP<const StridedMap> stridedRgFullMap = Teuchos::rcp_dynamic_cast<const StridedMap>(rangeMapExtractor->getFullMap());
331  Teuchos::RCP<const Map > fullRangeMap = Teuchos::null;
332  if(stridedRgFullMap != Teuchos::null) {
333  std::vector<size_t> stridedData = stridedRgFullMap->getStridingData();
334  fullRangeMap =
335  StridedMapFactory::Build(
336  bA->getRangeMap()->lib(),
338  fullRangeMapGIDs,
339  rangeIndexBase,
340  stridedData,
341  rebalancedComm, /*bA->getRangeMap()->getComm(),*/ //bA->getRangeMap()->getComm(),
342  stridedRgFullMap->getStridedBlockId(),
343  stridedRgFullMap->getOffset());
344  } else {
345  fullRangeMap =
346  MapFactory::Build(
347  bA->getRangeMap()->lib(),
349  fullRangeMapGIDs,
350  rangeIndexBase,
351  rebalancedComm /*bA->getRangeMap()->getComm()*/); //bA->getRangeMap()->getComm());
352  }
353  Teuchos::ArrayView<GO> fullDomainMapGIDs(fullDomainMapVector.size() ? &fullDomainMapVector[0] : 0,fullDomainMapVector.size());
354 
355  Teuchos::RCP<const StridedMap> stridedDoFullMap = Teuchos::rcp_dynamic_cast<const StridedMap>(domainMapExtractor->getFullMap());
356  Teuchos::RCP<const Map > fullDomainMap = Teuchos::null;
357  if(stridedDoFullMap != Teuchos::null) {
358  TEUCHOS_TEST_FOR_EXCEPTION(stridedDoFullMap==Teuchos::null, Exceptions::BadCast, "MueLu::RebalanceBlockedAc::Build: full map in domain map extractor has no striding information! error.");
359  std::vector<size_t> stridedData2 = stridedDoFullMap->getStridingData();
360  fullDomainMap =
361  StridedMapFactory::Build(
362  bA->getDomainMap()->lib(),
364  fullDomainMapGIDs,
365  domainIndexBase,
366  stridedData2,
367  rebalancedComm, /*bA->getDomainMap()->getComm(), *///bA->getDomainMap()->getComm(),
368  stridedDoFullMap->getStridedBlockId(),
369  stridedDoFullMap->getOffset());
370  } else {
371 
372  fullDomainMap =
373  MapFactory::Build(
374  bA->getDomainMap()->lib(),
376  fullDomainMapGIDs,
377  domainIndexBase,
378  rebalancedComm/*bA->getDomainMap()->getComm()*/); //bA->getDomainMap()->getComm());
379  }
380 
381  if(bRestrictComm) {
382  fullRangeMap->removeEmptyProcesses();
383  fullDomainMap->removeEmptyProcesses();
384  }
385 
386  // build map extractors
387  Teuchos::RCP<const MapExtractorClass> rebRangeMapExtractor = MapExtractorFactoryClass::Build(fullRangeMap, subBlockARangeMaps, bThyraRangeGIDs);
388  Teuchos::RCP<const MapExtractorClass> rebDomainMapExtractor = MapExtractorFactoryClass::Build(fullDomainMap, subBlockADomainMaps, bThyraDomainGIDs);
389 
390  TEUCHOS_TEST_FOR_EXCEPTION(rangeMapExtractor->NumMaps() != rebRangeMapExtractor->NumMaps(), Exceptions::BadCast, "MueLu::RebalanceBlockedAc::Build: Rebalanced RangeMapExtractor has " << rebRangeMapExtractor << " sub maps. Original RangeMapExtractor has " << rangeMapExtractor->NumMaps() << ". They must match!");
391  TEUCHOS_TEST_FOR_EXCEPTION(domainMapExtractor->NumMaps() != rebDomainMapExtractor->NumMaps(), Exceptions::BadCast, "MueLu::RebalanceBlockedAc::Build: Rebalanced DomainMapExtractor has " << rebDomainMapExtractor << " sub maps. Original DomainMapExtractor has " << domainMapExtractor->NumMaps() << ". They must match!");
392 
393  Teuchos::RCP<BlockedCrsMatrix> reb_bA = Teuchos::rcp(new BlockedCrsMatrix(rebRangeMapExtractor,rebDomainMapExtractor,10));
394  for(size_t i=0; i<bA->Rows(); i++) {
395  for(size_t j=0; j<bA->Cols(); j++) {
396  //Teuchos::RCP<const CrsMatrixWrap> crsOpij = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(subBlockRebA[i*bA->Cols() + j]);
397  reb_bA->setMatrix(i,j,subBlockRebA[i*bA->Cols() + j]);
398  }
399  }
400 
401  reb_bA->fillComplete();
402 
403  //reb_bA->describe(*out,Teuchos::VERB_EXTREME);
404  coarseLevel.Set("A", Teuchos::rcp_dynamic_cast<Matrix>(reb_bA), this);
405  // rebalance additional data:
406  // be aware, that we just call the rebalance factories without switching to local
407  // factory managers, i.e. the rebalance factories have to be defined with the appropriate
408  // factories by the user!
409  if (rebalanceFacts_.begin() != rebalanceFacts_.end()) {
410  SubFactoryMonitor m2(*this, "Rebalance additional data", coarseLevel);
411 
412  // call Build of all user-given transfer factories
413  for (std::vector<RCP<const FactoryBase> >::const_iterator it2 = rebalanceFacts_.begin(); it2 != rebalanceFacts_.end(); ++it2) {
414  GetOStream(Runtime0) << "RebalanceBlockedAc: call rebalance factory " << (*it2).get() << ": " << (*it2)->description() << std::endl;
415  (*it2)->CallBuild(coarseLevel);
416  }
417  }
418  } //Build()
419 
420  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
422 
423  /*TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::rcp_dynamic_cast<const TwoLevelFactoryBase>(factory) == Teuchos::null, Exceptions::BadCast,
424  "MueLu::RAPFactory::AddTransferFactory: Transfer factory is not derived from TwoLevelFactoryBase. "
425  "This is very strange. (Note: you can remove this exception if there's a good reason for)");
426  TEUCHOS_TEST_FOR_EXCEPTION(hasDeclaredInput_, Exceptions::RuntimeError, "MueLu::RAPFactory::AddTransferFactory: Factory is being added after we have already declared input");*/
427  rebalanceFacts_.push_back(factory);
428  } //AddRebalanceFactory()
429 
430 } //namespace MueLu
431 
432 #endif /* MUELU_REBALANCEBLOCKACFACTORY_DEF_HPP_ */
Exception indicating invalid cast attempted.
void AddFactoryManager(RCP< const FactoryManagerBase > FactManager)
Add a factory manager.
GlobalOrdinal GO
iterator begin() const
T & get(const std::string &name, T def_value)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
Timer to be used in factories. Similar to Monitor but with additional timers.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
#define SET_VALID_ENTRY(name)
Print additional debugging information.
One-liner description of what is happening.
Print statistics that do not involve significant additional computation.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void AddRebalanceFactory(const RCP< const FactoryBase > &factory)
Add rebalancing factory in the end of list of rebalancing factories in RebalanceAcFactory.
Class that holds all level-specific information.
Definition: MueLu_Level.hpp:99
Timer to be used in factories. Similar to SubMonitor but adds a timer level by level.
static std::string PrintMatrixInfo(const Matrix &A, const std::string &msgTag, RCP< const Teuchos::ParameterList > params=Teuchos::null)
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
iterator end() const
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
Exception throws to report errors in the internal logical of the program.
An exception safe way to call the method &#39;Level::SetFactoryManager()&#39;.
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()
bool is_null() const