MueLu  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node > Class Template Reference

Helper class for providing arbitrated communication across processors. More...

#include <MueLu_CoupledAggregationCommHelper_decl.hpp>

Inheritance diagram for MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >:
MueLu::BaseClass MueLu::VerboseObject MueLu::Describable Teuchos::VerboseObject< VerboseObject > Teuchos::Describable Teuchos::VerboseObjectBase Teuchos::LabeledObject

Public Member Functions

void ArbitrateAndCommunicate (Vector &weights, Aggregates &aggregates, const bool perturb) const
 This method assigns unknowns to aggregates. More...
 
void ArbitrateAndCommunicate (Vector &weight, LOVector &procWinner, LOMultiVector *companion, const bool perturb) const
 This class uses a weighted rendezvous algorithm to do a global reduction on a vector that may be based on a non unique map. More...
 
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 global id across many processors. More...
 
- Public Member Functions inherited from MueLu::BaseClass
virtual ~BaseClass ()
 Destructor. More...
 
- Public Member Functions inherited from MueLu::VerboseObject
VerbLevel GetVerbLevel () const
 Get the verbosity level. More...
 
void SetVerbLevel (const VerbLevel verbLevel)
 Set the verbosity level of this object. More...
 
int GetProcRankVerbose () const
 Get proc rank used for printing. Do not use this information for any other purpose. More...
 
int SetProcRankVerbose (int procRank) const
 Set proc rank used for printing. More...
 
bool IsPrint (MsgType type, int thisProcRankOnly=-1) const
 Find out whether we need to print out information for a specific message type. More...
 
Teuchos::FancyOStreamGetOStream (MsgType type, int thisProcRankOnly=0) const
 Get an output stream for outputting the input message type. More...
 
Teuchos::FancyOStreamGetBlackHole () const
 
 VerboseObject ()
 
virtual ~VerboseObject ()
 Destructor. More...
 
- Public Member Functions inherited from Teuchos::VerboseObject< VerboseObject >
 VerboseObject (const EVerbosityLevel verbLevel=VERB_DEFAULT, const RCP< FancyOStream > &oStream=Teuchos::null)
 
virtual const VerboseObjectsetVerbLevel (const EVerbosityLevel verbLevel) const
 
virtual const VerboseObjectsetOverridingVerbLevel (const EVerbosityLevel verbLevel) const
 
virtual EVerbosityLevel getVerbLevel () const
 
TEUCHOSPARAMETERLIST_LIB_DLL_EXPORT
RCP< const ParameterList
getValidVerboseObjectSublist ()
 
TEUCHOSPARAMETERLIST_LIB_DLL_EXPORT
void 
setupVerboseObjectSublist (ParameterList *paramList)
 
TEUCHOSPARAMETERLIST_LIB_DLL_EXPORT
void 
readVerboseObjectSublist (ParameterList *paramList, RCP< FancyOStream > *oStream, EVerbosityLevel *verbLevel)
 
void readVerboseObjectSublist (ParameterList *paramList, VerboseObject< ObjectType > *verboseObject)
 
- Public Member Functions inherited from Teuchos::VerboseObjectBase
virtual ~VerboseObjectBase ()
 
 VerboseObjectBase (const RCP< FancyOStream > &oStream=Teuchos::null)
 
virtual const VerboseObjectBasesetOStream (const RCP< FancyOStream > &oStream) const
 
virtual const VerboseObjectBasesetOverridingOStream (const RCP< FancyOStream > &oStream) const
 
virtual VerboseObjectBasesetLinePrefix (const std::string &linePrefix)
 
virtual RCP< FancyOStreamgetOStream () const
 
virtual RCP< FancyOStreamgetOverridingOStream () const
 
virtual std::string getLinePrefix () const
 
virtual OSTab getOSTab (const int tabs=1, const std::string &linePrefix="") const
 
- Public Member Functions inherited from MueLu::Describable
virtual ~Describable ()
 Destructor. More...
 
virtual std::string ShortClassName () const
 Return the class name of the object, without template parameters and without namespace. More...
 
virtual void describe (Teuchos::FancyOStream &out_arg, const VerbLevel verbLevel=Default) const
 
virtual std::string description () const
 Return a simple one-line description of this object. More...
 
void describe (Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
 Print the object with some verbosity level to an FancyOStream object. More...
 
- Public Member Functions inherited from Teuchos::Describable
void describe (std::ostream &out, const EVerbosityLevel verbLevel=verbLevel_default) const
 
 LabeledObject ()
 
virtual ~LabeledObject ()
 
virtual void setObjectLabel (const std::string &objectLabel)
 
virtual std::string getObjectLabel () const
 
DescribableStreamManipulatorState describe (const Describable &describable, const EVerbosityLevel verbLevel=Describable::verbLevel_default)
 
std::ostream & operator<< (std::ostream &os, const DescribableStreamManipulatorState &d)
 

Private Types

typedef DefaultScalar Scalar
 

Private Attributes

RCP< const Import > import_
 
RCP< const Import > winnerImport_
 
RCP< Import > pushWinners_
 
RCP< Vector > tempVec_
 
RCP< Vector > perturbWt_
 
RCP< Vector > postComm_
 
RCP< Vector > candidateWinners_
 
ArrayRCP< GOmyWinners_
 
int numMyWinners_
 
RCP< Map > winnerMap_
 
int numCalls_
 
int myPID_
 

Constructors/Destructors.

 CoupledAggregationCommHelper (const RCP< const Map > &uniqueMap, const RCP< const Map > &nonUniqueMap)
 Constructor. More...
 
 ~CoupledAggregationCommHelper ()
 Destructor. More...
 

Additional Inherited Members

- Static Public Member Functions inherited from MueLu::VerboseObject
static void SetMueLuOStream (const Teuchos::RCP< Teuchos::FancyOStream > &mueluOStream)
 
static void SetMueLuOFileStream (const std::string &filename)
 
static Teuchos::RCP
< Teuchos::FancyOStream
GetMueLuOStream ()
 
static void SetDefaultVerbLevel (const VerbLevel defaultVerbLevel)
 Set the default (global) verbosity level. More...
 
static VerbLevel GetDefaultVerbLevel ()
 Get the default (global) verbosity level. More...
 
- Static Public Member Functions inherited from Teuchos::VerboseObject< VerboseObject >
static void setDefaultVerbLevel (const EVerbosityLevel defaultVerbLevel)
 
static EVerbosityLevel getDefaultVerbLevel ()
 
- Static Public Member Functions inherited from Teuchos::VerboseObjectBase
static void setDefaultOStream (const RCP< FancyOStream > &defaultOStream)
 
static RCP< FancyOStreamgetDefaultOStream ()
 
- Static Public Attributes inherited from Teuchos::Describable
static const EVerbosityLevel verbLevel_default
 
- Protected Member Functions inherited from Teuchos::VerboseObject< VerboseObject >
void initializeVerboseObject (const EVerbosityLevel verbLevel=VERB_DEFAULT, const RCP< FancyOStream > &oStream=Teuchos::null)
 
- Protected Member Functions inherited from Teuchos::VerboseObjectBase
void initializeVerboseObjectBase (const RCP< FancyOStream > &oStream=Teuchos::null)
 
virtual void informUpdatedVerbosityState () const
 

Detailed Description

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
class MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >

Helper class for providing arbitrated communication across processors.

For more details, see the comments for the ArbitrateAndCommunicate methods.

Definition at line 73 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

Member Typedef Documentation

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
typedef DefaultScalar MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::Scalar
private

Definition at line 74 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

Constructor & Destructor Documentation

template<class LocalOrdinal , class GlobalOrdinal , class Node >
MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::CoupledAggregationCommHelper ( const RCP< const Map > &  uniqueMap,
const RCP< const Map > &  nonUniqueMap 
)

Constructor.

Definition at line 64 of file MueLu_CoupledAggregationCommHelper_def.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::~CoupledAggregationCommHelper ( )
inline

Destructor.

Definition at line 88 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

Member Function Documentation

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
void MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::ArbitrateAndCommunicate ( Vector &  weights,
Aggregates aggregates,
const bool  perturb 
) const
inline

This method assigns unknowns to aggregates.

Tie-breaking is possible is using random weights.

Parameters
[in]weightsvector of weights that help determine ownership.
[in,out]aggregatesaggregate data structure
[in]perturbflag indicating whether weights should be randomly perturbed for tie-breaking purposes.

Definition at line 103 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal , class GlobalOrdinal , class Node >
void MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::ArbitrateAndCommunicate ( Vector &  weight,
LOVector &  procWinner,
LOMultiVector *  companion,
const bool  perturb 
) const

This class uses a weighted rendezvous algorithm to do a global reduction on a vector that may be based on a non unique map.

A non-unique map is one that has at least one global ID that occurs on two or more processes. For each repeated ID \(i\), the algorithm finds the maximum value \(v[i]\) in the weight vector \(v\). This value is communicated to all processors that have \(i\) in their local map. More details are below.

For each GlobalId \(K\) associated with weight.getMap():

 -# Find the maximum absolute value of \form#4 across all
    processors and assign this to all local elements of weight[] (across
    processors) that are associated with \form#3.
 -# Set procWinner[] to the MyPid() that had the largest element.
    procWinner[] is still set if only one processor owns a GlobalId.

    The ONLY CASE when procWinner[i] is NOT set corresponds to when
    all local weights associated with a GlobalId are zero. This allows
    one to effectively skip the maximum/winner calculation for a subset
    of GlobalId's.  This might occur when a processor has already
    claimed ownership for a GlobalId and so all local copies have
    the same value. We want to skip the maximum calculation with
    tiebreaking to avoid another processor claiming ownership.

 -# Optionally, set companion[] (across all relevant processors) to the
    local companion value associated with the procWinner[] processor.

@param weight[in,out]
                        - On input, vector of NONNEGATIVE weights.
                        - On output, \form#5
                         where \form#6 is processor \form#7's value for GID \form#8.

@param procWinner[in,out]
                         - On input, allocated but contents ignored.
                         - On output, \form#9  such that

\(\mbox{weight}[k_{pj}] = \max(\mbox{weight}[k_{p1}],...,\mbox{weight}[k_{pn}])\), where \( \mbox{weight}[k_{pj}] \) is processor \(pj\)'s value for GID \(k\). NOTE: If all input \(\mbox{weight}[k_{pi}]\)'s are zero, then \(\mbox{procWinner}[k]\) is left untouched.

Parameters
companion[in,out]
  • On input, either NULL or allocated but contents ignored. If NULL, step 3 above is skipped.
  • On output, if not null, \(\mbox{companion}[k] \Leftarrow \mbox{companion}[k_j]\) where \(\mbox{companion}[k_j]\) lives on processor \(\mbox{procWinner}[k]\). and corresponds to the same GlobalId as \(k\). NOTE: If for a particular GlobalId, no processor has a value of procWinner that matches its MyPid, the corresponding companion is not altered.
 @param perturb[in]                  Optional arguments that is either true or
                                     false (default: true). weight is perturbed
                                     and the perturbed values are used in step 1)
                                     above. Returned values reflect the perturbed
                                     data. This option avoids having lots of
                                     tiebreaks where the large MyPid() always wins.

Definition at line 73 of file MueLu_CoupledAggregationCommHelper_def.hpp.

template<class LocalOrdinal , class GlobalOrdinal , class Node >
void MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::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 global id across many processors.

The source may not have the same value for all of these multiple copies, but on termination dest will have a unique value for each global id. When multiple copies exist in source, 'what' determines how they are combined to make a unique value in dest (see CombineMode).

Input:

Parameters
[in]sourceVector where multiple copies of some GlobalIds might exist and might have different values.
[in,out]destOn input, allocated but contents ignored. On output, contains redistributed data from source where 'what' determines how multiple copies of source values associated with the same GlobalId are combined into a unique value on all processors.
[in]whatDetermines how multiple copies of the same GlobalId are combined (see CombineMode).

Definition at line 475 of file MueLu_CoupledAggregationCommHelper_def.hpp.

Member Data Documentation

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
RCP<const Import> MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::import_
private

Definition at line 212 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
RCP<const Import> MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::winnerImport_
mutableprivate

Definition at line 213 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
RCP<Import> MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::pushWinners_
mutableprivate

Definition at line 214 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
RCP<Vector> MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::tempVec_
private

Definition at line 215 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
RCP<Vector> MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::perturbWt_
mutableprivate

Definition at line 216 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
RCP<Vector> MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::postComm_
mutableprivate

Definition at line 217 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
RCP<Vector> MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::candidateWinners_
mutableprivate

Definition at line 218 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
ArrayRCP<GO> MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::myWinners_
mutableprivate

Definition at line 219 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
int MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::numMyWinners_
mutableprivate

Definition at line 220 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
RCP<Map> MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::winnerMap_
mutableprivate

Definition at line 221 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
int MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::numCalls_
mutableprivate

Definition at line 222 of file MueLu_CoupledAggregationCommHelper_decl.hpp.

template<class LocalOrdinal = DefaultLocalOrdinal, class GlobalOrdinal = DefaultGlobalOrdinal, class Node = DefaultNode>
int MueLu::CoupledAggregationCommHelper< LocalOrdinal, GlobalOrdinal, Node >::myPID_
private

Definition at line 223 of file MueLu_CoupledAggregationCommHelper_decl.hpp.


The documentation for this class was generated from the following files: