MueLu  Version of the Day
MueLu_CoupledAggregationCommHelper_decl.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_COUPLEDAGGREGATIONCOMMHELPER_DECL_HPP
47 #define MUELU_COUPLEDAGGREGATIONCOMMHELPER_DECL_HPP
48 
49 #include <Xpetra_Import_fwd.hpp>
51 #include <Xpetra_Vector_fwd.hpp>
53 
54 #include "MueLu_ConfigDefs.hpp"
55 #include "MueLu_BaseClass.hpp"
57 
58 #include "MueLu_Aggregates.hpp"
59 
60 namespace MueLu {
61 
69  template <class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
71 
72  typedef double Scalar; // Scalar type only used for weight: always a double.
73 #undef MUELU_COUPLEDAGGREGATIONCOMMHELPER_SHORT
74 #include "MueLu_UseShortNames.hpp"
75 
76  public:
77 
79 
80 
82  CoupledAggregationCommHelper(const RCP<const Map> & uniqueMap, const RCP<const Map> & nonUniqueMap);
83 
86 
88 
100  void ArbitrateAndCommunicate(Vector &weights, Aggregates &aggregates, const bool perturb) const {
101  ArbitrateAndCommunicate(weights, *aggregates.GetProcWinner(), &*aggregates.GetVertex2AggId(), perturb);
102  }
103 
162  /*
163  Output:
164  @param weight \f$ weight[k] \Leftarrow \max(weight[k_1],\dots,weight[k_n]) \f$
165  where \f$ weight[k_j] \f$ live on different processors
166  but have the same GlobalId as weight[k] on this processor.
167 
168  @param procWinner procWinner[k] <-- MyPid associated with the
169  kj yielding the max in
170  Max(weight[k1],...,weight[kn]) .
171  See weight Output comments.
172  NOTE: If all input weight[kj]'s are zero,
173  then procWinner[k] is left untouched.
174 
175  @param companion If not null,
176  companion[k] <-- companion[kj] where
177  companion[kj] lives on processor procWinner[k].
178  and corresponds to the same GlobalId as k.
179  NOTE: If for a particlar GlobalId, no processor
180  has a value of procWinner that matches
181  its MyPid, the corresponding companion
182  is not altered.
183  */
184  void ArbitrateAndCommunicate(Vector &weight, LOVector &procWinner, LOMultiVector *companion, const bool perturb) const; //ArbitrateAndCommunicate(Vector&, LOVector &, LOMultiVector *, const bool) const
185 
206  void NonUnique2NonUnique(const Vector &source, Vector &dest, const Xpetra::CombineMode what) const;
207 
208  private:
210  mutable RCP<const Import> winnerImport_; //FIXME get rid of "mutable"
211  mutable RCP<Import> pushWinners_; //FIXME get rid of mutable
217  mutable int numMyWinners_;
219  mutable int numCalls_;
220  int myPID_;
221 
222  // uniqueMap A subset of weight.getMap() where each GlobalId
223  // has only one unique copy on one processor.
224  // Normally, weight.getMap() would have both locals
225  // and ghost elements while uniqueMap would just
226  // have the locals. It should be possible to
227  // remove this or make it an optional argument
228  // and use some existing Epetra/Tpetra capability to
229  // make a uniqueMap.
230  //
231  // import_ This corresponds precisely to
232  // Import import_(
233  // weight.getMap(), uniqueMap);
234  // This could also be eliminated and created
235  // here, but for efficiency user's must pass in.
236  //
237  };
238 
239 
240 }
241 
242 //JG:
243 // - procWinner is an array of proc ID -> LocalOrdinal
244 // - companion == aggregates.GetVertex2AggId() == local aggregate ID -> LocalOrdinal
245 
246 #define MUELU_COUPLEDAGGREGATIONCOMMHELPER_SHORT
247 #endif // MUELU_COUPLEDAGGREGATIONCOMMHELPER_DECL_HPP
Container class for aggregation information.
Namespace for MueLu classes and methods.
void NonUnique2NonUnique(const Vector &source, Vector &dest, const Xpetra::CombineMode what) const
Redistribute data in source to dest where both source and dest might have multiple copies of the same...
Helper class for providing arbitrated communication across processors.
Base class for MueLu classes.
CoupledAggregationCommHelper(const RCP< const Map > &uniqueMap, const RCP< const Map > &nonUniqueMap)
Constructor.
const RCP< LOMultiVector > & GetVertex2AggId() const
Returns constant vector that maps local node IDs to local aggregates IDs.
const RCP< LOVector > & GetProcWinner() const
Returns constant vector that maps local node IDs to owning processor IDs.
void ArbitrateAndCommunicate(Vector &weights, Aggregates &aggregates, const bool perturb) const
This method assigns unknowns to aggregates.