MueLu Version of the Day
Loading...
Searching...
No Matches
MueLu_Aggregates_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_AGGREGATES_DEF_HPP
47#define MUELU_AGGREGATES_DEF_HPP
48
49#include <Xpetra_Map.hpp>
50#include <Xpetra_Vector.hpp>
51#include <Xpetra_MultiVector.hpp>
52#include <Xpetra_BlockedMultiVector.hpp>
53#include <Xpetra_BlockedVector.hpp>
54#include <Xpetra_VectorFactory.hpp>
55#include <Xpetra_MultiVectorFactory.hpp>
56
58#include "MueLu_Graph.hpp"
59#include "MueLu_Utilities_decl.hpp" // MueLu_sumAll
60
61namespace MueLu {
62
64 template <class LocalOrdinal, class GlobalOrdinal, class Node>
66 nAggregates_ = 0;
67 nGlobalAggregates_ = 0;
68
69 vertex2AggId_ = LOMultiVectorFactory::Build(graph.GetImportMap(), 1);
70 vertex2AggId_->putScalar(MUELU_UNAGGREGATED);
71
72 procWinner_ = LOVectorFactory::Build(graph.GetImportMap());
73 procWinner_->putScalar(MUELU_UNASSIGNED);
74
75 isRoot_ = Teuchos::ArrayRCP<bool>(graph.GetImportMap()->getLocalNumElements(), false);
76
77 // slow but safe, force TentativePFactory to build column map for P itself
78 aggregatesIncludeGhosts_ = true;
79 }
80
82 template <class LocalOrdinal, class GlobalOrdinal, class Node>
84 nAggregates_ = 0;
85 nGlobalAggregates_ = 0;
86
87 vertex2AggId_ = LOMultiVectorFactory::Build(map, 1);
88 vertex2AggId_->putScalar(MUELU_UNAGGREGATED);
89
90 procWinner_ = LOVectorFactory::Build(map);
91 procWinner_->putScalar(MUELU_UNASSIGNED);
92
93 isRoot_ = Teuchos::ArrayRCP<bool>(map->getLocalNumElements(), false);
94
95 // slow but safe, force TentativePFactory to build column map for P itself
96 aggregatesIncludeGhosts_ = true;
97 }
98
100 template <class LocalOrdinal, class GlobalOrdinal, class Node>
101 Teuchos::ArrayRCP<LocalOrdinal> Aggregates<LocalOrdinal, GlobalOrdinal, Node>::ComputeAggregateSizes(bool forceRecompute) const {
102
103 if (aggregateSizes_ != Teuchos::null && !forceRecompute) {
104
105 return aggregateSizes_;
106
107 } else {
108
109 //invalidate previous sizes.
110 aggregateSizes_ = Teuchos::null;
111
112 Teuchos::ArrayRCP<LO> aggregateSizes;
113 aggregateSizes = Teuchos::ArrayRCP<LO>(nAggregates_,0);
114 int myPid = vertex2AggId_->getMap()->getComm()->getRank();
115 Teuchos::ArrayRCP<LO> procWinner = procWinner_->getDataNonConst(0);
116 Teuchos::ArrayRCP<LO> vertex2AggId = vertex2AggId_->getDataNonConst(0);
117 LO size = procWinner.size();
118
119 //for (LO i = 0; i < nAggregates_; ++i) aggregateSizes[i] = 0;
120 for (LO k = 0; k < size; ++k ) {
121 if (procWinner[k] == myPid) aggregateSizes[vertex2AggId[k]]++;
122 }
123
124 aggregateSizes_ = aggregateSizes;
125
126 return aggregateSizes;
127 }
128
129 } //ComputeAggSizesNodes
130
131 template <class LocalOrdinal, class GlobalOrdinal, class Node>
133 std::ostringstream out;
134 out << BaseClass::description();
135 if (nGlobalAggregates_ == -1) out << "{nGlobalAggregates = not computed}";
136 else out << "{nGlobalAggregates = " << nGlobalAggregates_ << "}";
137 return out.str();
138 }
139
140 template <class LocalOrdinal, class GlobalOrdinal, class Node>
141 void Aggregates<LocalOrdinal, GlobalOrdinal, Node>::print(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const {
143
144
145 if (verbLevel & Statistics1) {
146 if (nGlobalAggregates_ == -1) out0 << "Global number of aggregates: not computed " << std::endl;
147 else out0 << "Global number of aggregates: " << nGlobalAggregates_ << std::endl;
148 }
149
150 if(verbLevel == Teuchos::VERB_EXTREME) {
151 for(size_t j=0; j <vertex2AggId_->getNumVectors(); j++) {
152 auto data = vertex2AggId_->getData(j);
153 for(size_t i=0; i<vertex2AggId_->getLocalLength(); i++) {
154 out<<i<<" : "<< data[i] <<std::endl;
155 }
156 out<<std::endl;
157 }
158 }
159
160 }
161
162 template <class LocalOrdinal, class GlobalOrdinal, class Node>
164 if (nGlobalAggregates_ != -1) {
165 LO nAggregates = GetNumAggregates();
166 GO nGlobalAggregates;
167 MueLu_sumAll(vertex2AggId_->getMap()->getComm(), (GO)nAggregates, nGlobalAggregates);
168 SetNumGlobalAggregates(nGlobalAggregates);
169 }
170 return nGlobalAggregates_;
171 }
172
173 template <class LocalOrdinal, class GlobalOrdinal, class Node>
174 const RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal, Node> > Aggregates<LocalOrdinal, GlobalOrdinal, Node>::GetMap() const {
175 return vertex2AggId_->getMap();
176 }
177
178 template <class LocalOrdinal, class GlobalOrdinal, class Node>
179 void Aggregates<LocalOrdinal, GlobalOrdinal, Node>::ComputeNodesInAggregate(Array<LO> & aggPtr, Array<LO> & aggNodes, Array<LO> & unaggregated) const {
180 LO numAggs = GetNumAggregates();
181 LO numNodes = vertex2AggId_->getLocalLength();
182 Teuchos::ArrayRCP<const LO> vertex2AggId = vertex2AggId_->getData(0);
183 Teuchos::ArrayRCP<LO> aggSizes = ComputeAggregateSizes(true);
184 LO INVALID = Teuchos::OrdinalTraits<LO>::invalid();
185
186 aggPtr.resize(numAggs+1);
187 Array<LO> aggCurr(numAggs+1);
188 aggNodes.resize(numNodes);
189 unaggregated.resize(numNodes);
190
191 LO currNumUnaggregated=0;
192
193 // Construct the "rowptr" and the counter
194 aggPtr[0] = 0;
195 for(LO i=0; i<numAggs; i++) {
196 aggPtr[i+1] = aggSizes[i] + aggPtr[i];
197 aggCurr[i] = aggPtr[i];
198 }
199
200 // Stick the nodes in each aggregate's spot
201 for(LO i=0; i<numNodes; i++) {
202 LO aggregate = vertex2AggId[i];
203 if(aggregate !=INVALID) {
204 aggNodes[aggCurr[aggregate]] = i;
205 aggCurr[aggregate]++;
206 }
207 else {
208 unaggregated[currNumUnaggregated] = i;
209 currNumUnaggregated++;
210 }
211 }
212 unaggregated.resize(currNumUnaggregated);
213
214 }
215
216
217} //namespace MueLu
218
219#endif // MUELU_AGGREGATES_DEF_HPP
#define MUELU_UNAGGREGATED
#define MUELU_UNASSIGNED
#define MUELU_DESCRIBE
Helper macro for implementing Describable::describe() for BaseClass objects.
#define MueLu_sumAll(rcpComm, in, out)
MueLu::DefaultGlobalOrdinal GlobalOrdinal
Teuchos::ArrayRCP< LO > ComputeAggregateSizes(bool forceRecompute=false) const
Compute sizes of aggregates.
Aggregates(const GraphBase &graph)
Standard constructor for Aggregates structure.
void print(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
const RCP< const Map > GetMap() const
returns (overlapping) map of aggregate/node distribution
GO GetNumGlobalAggregatesComputeIfNeeded()
Get global number of aggregates.
std::string description() const
Return a simple one-line description of this object.
void ComputeNodesInAggregate(Array< LO > &aggPtr, Array< LO > &aggNodes, Array< LO > &unaggregated) const
Generates a compressed list of nodes in each aggregate, where the entries in aggNodes[aggPtr[i]] up t...
virtual std::string description() const
Return a simple one-line description of this object.
MueLu representation of a graph.
virtual const RCP< const Map > GetImportMap() const =0
Namespace for MueLu classes and methods.
@ Statistics1
Print more statistics.