Zoltan2
Loading...
Searching...
No Matches
Zoltan2_GraphMetricsUtility.hpp
Go to the documentation of this file.
1// @HEADER
2//
3// ***********************************************************************
4//
5// Zoltan2: A package of combinatorial algorithms for scientific computing
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 Karen Devine (kddevin@sandia.gov)
39// Erik Boman (egboman@sandia.gov)
40// Siva Rajamanickam (srajama@sandia.gov)
41//
42// ***********************************************************************
43//
44// @HEADER
45
49#ifndef ZOLTAN2_GRAPHICMETRICVALUESUTILITY_HPP
50#define ZOLTAN2_GRAPHICMETRICVALUESUTILITY_HPP
51
55#include <zoltan_dd.h>
56#include <Zoltan2_TPLTraits.hpp>
57#include <map>
58#include <vector>
59
60namespace Zoltan2 {
61
90template <typename Adapter,
91 typename MachineRep = // Default MachineRep type
92 MachineRepresentation<typename Adapter::scalar_t,
93 typename Adapter::part_t> >
95 const RCP<const Environment> &env,
96 const RCP<const Comm<int> > &comm,
97 const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graph,
98 const ArrayView<const typename Adapter::part_t> &parts,
99 typename Adapter::part_t &numParts,
100 ArrayRCP<RCP<BaseClassMetrics<typename Adapter::scalar_t> > > &metrics,
101 ArrayRCP<typename Adapter::scalar_t> &globalSums,
102 bool bMessages = true,
103 const RCP <const MachineRep> machine = Teuchos::null) {
104
105 env->timerStart(MACRO_TIMERS, "globalWeightedByPart");
106
107 // Note we used to have with hops as a separate method but decided to combine
108 // both into this single method. machine is an optional parameter to choose
109 // between the two methods.
110 bool bHops = (machine != Teuchos::null);
111
112 env->debug(DETAILED_STATUS, "Entering globalWeightedByPart");
113
115 // Initialize return values
116
117 typedef typename Adapter::lno_t lno_t;
118 typedef typename Adapter::gno_t gno_t;
119 typedef typename Adapter::offset_t offset_t;
120 typedef typename Adapter::scalar_t scalar_t;
121 typedef typename Adapter::part_t part_t;
122
124 t_input_t;
125
126 lno_t localNumVertices = graph->getLocalNumVertices();
127 lno_t localNumEdges = graph->getLocalNumEdges();
128
129 ArrayView<const gno_t> Ids;
130 ArrayView<t_input_t> v_wghts;
131 graph->getVertexList(Ids, v_wghts);
132
133 typedef GraphMetrics<scalar_t> gm_t;
134
135 // get the edge ids, and weights
136 ArrayView<const gno_t> edgeIds;
137 ArrayView<const offset_t> offsets;
138 ArrayView<t_input_t> e_wgts;
139 graph->getEdgeList(edgeIds, offsets, e_wgts);
140
141
142 std::vector <scalar_t> edge_weights;
143 int numWeightPerEdge = graph->getNumWeightsPerEdge();
144
145 int numMetrics = bHops ?
146 4 : // "edge cuts", messages, hops, weighted hops
147 2; // "edge cuts", messages
148
149 if (numWeightPerEdge) numMetrics += bHops ?
150 numWeightPerEdge * 2 : // "weight n", weighted hops per weight n
151 numWeightPerEdge; // "weight n"
152
153 // add some more metrics to the array
154 auto next = metrics.size(); // where we begin filling
155 for (auto n = 0; n < numMetrics; ++n) {
156 addNewMetric<gm_t, scalar_t>(env, metrics);
157 }
158
159 std::vector <part_t> e_parts (localNumEdges);
160
161 std::vector<part_t> local_parts;
162
163#ifdef HAVE_ZOLTAN2_MPI
164 if (comm->getSize() > 1) {
165 const bool bUseLocalIDs = false; // Local IDs not needed
167 int debug_level = 0;
168 directory_t directory(comm, bUseLocalIDs, debug_level);
169
170 if (localNumVertices)
171 directory.update(localNumVertices, &Ids[0], NULL, &parts[0],
172 NULL, directory_t::Update_Mode::Replace);
173 else
174 directory.update(localNumVertices, NULL, NULL, NULL,
175 NULL, directory_t::Update_Mode::Replace);
176
177 if (localNumEdges)
178 directory.find(localNumEdges, &edgeIds[0], NULL, &e_parts[0],
179 NULL, NULL, false);
180 else
181 directory.find(localNumEdges, NULL, NULL, NULL,
182 NULL, NULL, false);
183 } else
184#endif
185 {
186 std::map<gno_t,lno_t> global_id_to_local_index;
187
188 // else everything is local.
189 // we need a globalid to local index conversion.
190 // this does not exists till this point, so we need to create one.
191
192 for (lno_t i = 0; i < localNumVertices; ++i){
193 //at the local index i, we have the global index Ids[i].
194 //so write i, to Ids[i] index of the vector.
195 global_id_to_local_index[Ids[i]] = i;
196 }
197
198 for (lno_t i = 0; i < localNumEdges; ++i){
199 gno_t ei = edgeIds[i];
200 //ei is the global index of the neighbor one.
201 part_t p = parts[global_id_to_local_index[ei]];
202 e_parts[i] = p;
203 }
204 }
205
206 RCP<const Teuchos::Comm<int> > tcomm = comm;
207
208 env->timerStart(MACRO_TIMERS, "Communication Graph Create");
209 {
210 const bool bUseLocalIDs = false; // Local IDs not needed
211 int debug_level = 0;
212 // Create a directory indexed by part_t with values t_scalar_t for weight sums
213
214 // this struct is the user data type for a part
215 // each part will have a std::vector<part_info> for its user data
216 // representing the list of all neighbors and a weight for each.
217 struct part_info {
218 part_info() : sum_weights(0) {
219 }
220
221 // operator +=
222 // this allows the directory to know how to assemble two structs
223 // which return true for ==.
224 // TODO: Decide if we want directory to work like this for AggregateAdd
225 const part_info & operator+=(const part_info & src) {
226 sum_weights += src.sum_weights;
227 return *this; // return old value
228 }
229
230 // operator>
231 // TODO: Decide if we want directory to work like this for AggregateAdd
232 bool operator>(const part_info & src) {
233 // Note: Currently this doesn't actually do anything except allow this
234 // struct to compile. Aggregate mode used operator> to preserve ordering
235 // and therefore a custom struct must currently define it. However in
236 // this test we will be using AggregateAdd mode which doesn't actually
237 // use this operator> . However if we change the test so we require the
238 // input data to already be ordered by target_part, then we could change
239 // the directory implementation so AggregateAdd and Aggregate are almost
240 // identical. The only difference would be that Aggregate mode would
241 // check operator== and if true, throw away the duplicate, while
242 // AggregateAdd mode would check operator== and if true, call operator+=
243 // to combine the values, in this case summing sum_weights.
244 return (target_part > src.target_part);
245 }
246
247 // operator==
248 // this allows the directory to know that two structs should be
249 // aggregated into one using the operator +=.
250 // TODO: Decide if we want directory to work like this for AggregateAdd
251 // This works but seems fussy/complicated - to discuss. I'm not yet sure
252 // how to best integrate this so we can aggregate both simple types where
253 // we just keep unique elements and more complex structs with a 'rule'
254 // for combining them.
255 bool operator==(const part_info & src) {
256 // if target_part is the same then the values for sum_weights will
257 // be summed.
258 return (target_part == src.target_part);
259 }
260
261 part_t target_part; // the part this part_info refers to
262 scalar_t sum_weights; // the sum of weights
263 };
264
265 // get the vertices in each part in my part.
266 std::vector <lno_t> part_begins(numParts, -1);
267 std::vector <lno_t> part_nexts(localNumVertices, -1);
268
269 // cluster vertices according to their parts.
270 // create local part graph.
271 for (lno_t i = 0; i < localNumVertices; ++i){
272 part_t ap = parts[i];
273 part_nexts[i] = part_begins[ap];
274 part_begins[ap] = i;
275 }
276
277 for (int weight_index = -1;
278 weight_index < numWeightPerEdge; ++weight_index) {
279
280 std::vector<part_t> part_data(numParts); // will resize to lower as needed
281 std::vector<std::vector<part_info>> user_data(numParts); // also to resize
282 int count_total_entries = 0;
283
284 std::vector <part_t> part_neighbors(numParts);
285 std::vector <scalar_t> part_neighbor_weights(numParts, 0);
286 std::vector <scalar_t> part_neighbor_weights_ordered(numParts);
287
288 // coarsen for all vertices in my part in order with parts.
289 for (part_t i = 0; i < numParts; ++i) {
290 part_t num_neighbor_parts = 0;
291 lno_t v = part_begins[i];
292 // get part i, and first vertex in this part v.
293 while (v != -1){
294 // now get the neightbors of v.
295 for (offset_t j = offsets[v]; j < offsets[v+1]; ++j){
296
297 // get the part of the second vertex.
298 part_t ep = e_parts[j];
299
300 // TODO: Can we skip condition (i==ep)
301 // The self reference set is going to be excluded later anyways
302 // so we could make this more efficient.
303 scalar_t ew = 1;
304 if (weight_index > -1){
305 ew = e_wgts[weight_index][j];
306 }
307 // add it to my local part neighbors for part i.
308 if (part_neighbor_weights[ep] < 0.00001){
309 part_neighbors[num_neighbor_parts++] = ep;
310 }
311 part_neighbor_weights[ep] += ew;
312 }
313 v = part_nexts[v];
314 }
315
316 // now get the part list.
317 for (lno_t j = 0; j < num_neighbor_parts; ++j) {
318 part_t neighbor_part = part_neighbors[j];
319 part_neighbor_weights_ordered[j] =
320 part_neighbor_weights[neighbor_part];
321 part_neighbor_weights[neighbor_part] = 0;
322 }
323
324 if (num_neighbor_parts > 0) {
325 // for the new directory note a difference in the logic flow
326 // originally we have CrsMatrix which could collect these values
327 // as we built each row. For the directory it's probably better to
328 // have update called just once so we collect the values and then
329 // do all of the update at the end.
330 part_data[count_total_entries] = i; // set up for directory
331 std::vector<part_info> & add_user_data =
332 user_data[count_total_entries];
333 ++count_total_entries;
334
335 add_user_data.resize(num_neighbor_parts);
336
337 for(int n = 0; n < num_neighbor_parts; ++n) {
338 part_info & add_data = add_user_data[n];
339 add_data.target_part = part_neighbors[n];
340 add_data.sum_weights = part_neighbor_weights_ordered[n];
341 }
342 }
343 }
344
345 scalar_t max_edge_cut = 0;
346 scalar_t total_edge_cut = 0;
347 part_t max_message = 0;
348 part_t total_message = 0;
349
350 part_t total_hop_count = 0;
351 scalar_t total_weighted_hop_count = 0;
352 part_t max_hop_count = 0;
353 scalar_t max_weighted_hop_count = 0;
354
355 // for serial or comm size 1 we need to fill this from local data
356 // TODO: Maybe remove all special casing for serial and make this pipeline
357 // uniform always
358 if(local_parts.size() == 0) {
359 local_parts.resize(numParts);
360 for(size_t n = 0; n < local_parts.size(); ++n) {
361 local_parts[n] = n;
362 }
363 }
364
365 std::vector<std::vector<part_info>> find_user_data;
366
367 // directory does not yet support SerialComm because it still has older
368 // MPI calls which need to be refactored to Teuchos::comm format. To
369 // work around this issue skip the directory calls and just set the
370 // find data equal to the input update data. This works because above
371 // logic has already summed the weights per process so in the SerialComm
372 // case, there won't be duplicates.
373 bool bSerialComm =
374 (dynamic_cast<const Teuchos::SerialComm<int>*>(&(*comm)) != NULL);
375
376 if(!bSerialComm) {
378 directory_t;
379 directory_t directory(comm, bUseLocalIDs, debug_level);
380
381 if(count_total_entries) {
382 // update
383 directory.update(count_total_entries, &part_data[0],
384 NULL, &user_data[0],
385 NULL, directory_t::Update_Mode::AggregateAdd);
386 }
387 else {
388 directory.update(count_total_entries, NULL, NULL, NULL,
389 NULL, directory_t::Update_Mode::AggregateAdd);
390 }
391
392 // get my local_parts (parts managed on this directory)
393 directory.get_locally_managed_gids(local_parts);
394
395 // set up find_user_data to have the right size
396 find_user_data.resize(local_parts.size());
397
398 // find
399 directory.find(local_parts.size(), &local_parts[0], NULL,
400 &find_user_data[0], NULL, NULL, false);
401 }
402 else {
403 find_user_data = user_data;
404 }
405
406 for(size_t n = 0; n < local_parts.size(); ++n) {
407 scalar_t part_edge_cut = 0;
408 part_t part_messages = 0;
409 const std::vector<part_info> & data = find_user_data[n];
410 for(size_t q = 0; q < data.size(); ++q) {
411 const part_t local_part = local_parts[n];
412 const part_info & info = data[q];
413 if (info.target_part != local_part) {
414 part_edge_cut += info.sum_weights;
415 part_messages += 1;
416
417 if(bHops) {
418 typename MachineRep::machine_pcoord_t hop_count = 0;
419 machine->getHopCount(local_part, info.target_part, hop_count);
420
421 part_t hop_counts = hop_count;
422 scalar_t weighted_hop_counts = hop_count * info.sum_weights;
423
424 total_hop_count += hop_counts;
425 total_weighted_hop_count += weighted_hop_counts;
426
427 if (hop_counts > max_hop_count ){
428 max_hop_count = hop_counts;
429 }
430 if (weighted_hop_counts > max_weighted_hop_count ){
431 max_weighted_hop_count = weighted_hop_counts;
432 }
433 }
434 }
435 }
436
437 if(part_edge_cut > max_edge_cut) {
438 max_edge_cut = part_edge_cut;
439 }
440 total_edge_cut += part_edge_cut;
441
442 if (part_messages > max_message){
443 max_message = part_messages;
444 }
445 total_message += part_messages;
446 }
447
448
449 scalar_t g_max_edge_cut = 0;
450 scalar_t g_total_edge_cut = 0;
451 part_t g_max_message = 0;
452 part_t g_total_message = 0;
453
454 part_t g_total_hop_count = 0;
455 scalar_t g_total_weighted_hop_count = 0;
456 part_t g_max_hop_count = 0;
457 scalar_t g_max_weighted_hop_count = 0;
458
459 try{
460 Teuchos::reduceAll<int, scalar_t>(*comm, Teuchos::REDUCE_MAX, 1,
461 &max_edge_cut, &g_max_edge_cut);
462 Teuchos::reduceAll<int, part_t>(*comm, Teuchos::REDUCE_MAX, 1,
463 &max_message, &g_max_message);
464
465 Teuchos::reduceAll<int, scalar_t>(*comm, Teuchos::REDUCE_SUM, 1,
466 &total_edge_cut, &g_total_edge_cut);
467 Teuchos::reduceAll<int, part_t>(*comm, Teuchos::REDUCE_SUM, 1,
468 &total_message, &g_total_message);
469
470 if(bHops) {
471 Teuchos::reduceAll<int, part_t>(*comm, Teuchos::REDUCE_MAX, 1,
472 &max_hop_count, &g_max_hop_count);
473 Teuchos::reduceAll<int, scalar_t>(*comm, Teuchos::REDUCE_MAX, 1,
474 &max_weighted_hop_count,
475 &g_max_weighted_hop_count);
476
477 Teuchos::reduceAll<int, part_t>(*comm, Teuchos::REDUCE_SUM, 1,
478 &total_hop_count, &g_total_hop_count);
479 Teuchos::reduceAll<int, scalar_t>(*comm, Teuchos::REDUCE_SUM, 1,
480 &total_weighted_hop_count,
481 &g_total_weighted_hop_count);
482 }
483 }
485
486 if (weight_index == -1){
487 metrics[next]->setName("edge cuts");
488 }
489 else {
490 std::ostringstream oss;
491 oss << "weight " << weight_index;
492 metrics[next]->setName( oss.str());
493 }
494 metrics[next]->setMetricValue("global maximum", g_max_edge_cut);
495 metrics[next]->setMetricValue("global sum", g_total_edge_cut);
496 next++;
497 if (weight_index == -1){
498 metrics[next]->setName("message");
499 metrics[next]->setMetricValue("global maximum", g_max_message);
500 metrics[next]->setMetricValue("global sum", g_total_message);
501 next++;
502 }
503
504 if(bHops) {
505 if (weight_index == -1){
506 metrics[next]->setName("hops (No Weight)");
507 metrics[next]->setMetricValue("global maximum", g_max_hop_count);
508 metrics[next]->setMetricValue("global sum", g_total_hop_count);
509 next++;
510 }
511
512 std::ostringstream oss;
513 oss << "weighted hops" << weight_index;
514 metrics[next]->setName( oss.str());
515 metrics[next]->
516 setMetricValue("global maximum", g_max_weighted_hop_count);
517 metrics[next]->
518 setMetricValue("global sum", g_total_weighted_hop_count);
519 next++;
520 }
521 }
522 }
523
524 env->timerStop(MACRO_TIMERS, "globalWeightedByPart");
525
526 env->debug(DETAILED_STATUS, "Exiting globalWeightedByPart");
527}
528
531template <typename scalar_t, typename part_t>
532void printGraphMetricsHeader(std::ostream &os,
533 part_t targetNumParts,
534 part_t numParts) {
535
536 os << "Graph Metrics: (" << numParts << " parts)";
537 os << std::endl;
538 if (targetNumParts != numParts) {
539 os << "Target number of parts is: " << targetNumParts << std::endl;
540 }
542}
543
546template <typename scalar_t, typename part_t>
547void printGraphMetrics(std::ostream &os,
548 part_t targetNumParts,
549 part_t numParts,
550 const ArrayView<RCP<BaseClassMetrics<scalar_t> > >
551 &infoList) {
552
553 printGraphMetricsHeader<scalar_t, part_t>(os, targetNumParts, numParts);
554 for (int i=0; i < infoList.size(); i++) {
555 if (infoList[i]->getName() != METRICS_UNSET_STRING) {
556 infoList[i]->printLine(os);
557 }
558 }
559 os << std::endl;
560}
561
564template <typename scalar_t, typename part_t>
565void printGraphMetrics(std::ostream &os,
566 part_t targetNumParts,
567 part_t numParts,
568 RCP<BaseClassMetrics<scalar_t>> metricValue) {
569
570 printGraphMetricsHeader<scalar_t, part_t>(os, targetNumParts, numParts);
571 metricValue->printLine(os);
572}
573
574} //namespace Zoltan2
575
576#endif
#define METRICS_UNSET_STRING
#define Z2_THROW_OUTSIDE_ERROR(env)
Throw an error returned from outside the Zoltan2 library.
Traits class to handle conversions between gno_t/lno_t and TPL data types (e.g., ParMETIS's idx_t,...
static void printHeader(std::ostream &os)
Print a standard header.
GraphModel defines the interface required for graph models.
size_t getLocalNumVertices() const
Returns the number vertices on this process.
map_t::local_ordinal_type lno_t
Definition: mapRemotes.cpp:17
map_t::global_ordinal_type gno_t
Definition: mapRemotes.cpp:18
Created by mbenlioglu on Aug 31, 2020.
void printGraphMetrics(std::ostream &os, part_t targetNumParts, part_t numParts, const ArrayView< RCP< BaseClassMetrics< scalar_t > > > &infoList)
Print out list of graph metrics.
@ MACRO_TIMERS
Time an algorithm (or other entity) as a whole.
void printGraphMetricsHeader(std::ostream &os, part_t targetNumParts, part_t numParts)
Print out header info for graph metrics.
@ DETAILED_STATUS
sub-steps, each method's entry and exit
void globalWeightedByPart(const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, const RCP< const GraphModel< typename Adapter::base_adapter_t > > &graph, const ArrayView< const typename Adapter::part_t > &parts, typename Adapter::part_t &numParts, ArrayRCP< RCP< BaseClassMetrics< typename Adapter::scalar_t > > > &metrics, ArrayRCP< typename Adapter::scalar_t > &globalSums, bool bMessages=true, const RCP< const MachineRep > machine=Teuchos::null)
Given the local partitioning, compute the global weighted cuts in each part.
SparseMatrixAdapter_t::part_t part_t