Zoltan2
Loading...
Searching...
No Matches
Metric.cpp
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//
46// Test the following:
47// EvaluatePartition class
48// MetricValues class
49// Metric related namespace methods
50
51
56#include <stdlib.h>
57#include <vector>
58
59
60using Teuchos::ArrayRCP;
61using Teuchos::Array;
62using Teuchos::RCP;
63using Teuchos::rcp;
64using Teuchos::arcp;
65
66using namespace std;
67using std::endl;
68using std::cout;
69
70template<class idInput_t>
71void doTest(RCP<const Comm<int> > comm, int numLocalObj,
72 int nWeights, int numLocalParts, bool givePartSizes, bool useDegreeAsWeight=false);
73
75
76// for testing basic input adapter
78
79// for testing graph adapter
80typedef Tpetra::CrsGraph<zlno_t, zgno_t, znode_t> tcrsGraph_t;
82
83// creates this so we can run the test suite over BasicIdentifierAdapter
84// and XpetraCrsGraphAdapter
85template<class idInput_t> void runTestSuite(RCP<const Comm<int> > comm, bool bCanTestDegreeAsWeights) {
86 doTest<idInput_t>(comm, 10, 0, -1, false);
87 doTest<idInput_t>(comm, 10, 0, 1, false);
88 doTest<idInput_t>(comm, 10, 0, 1, true);
89 doTest<idInput_t>(comm, 10, 1, 1, false);
90 doTest<idInput_t>(comm, 10, 1, 1, true);
91 doTest<idInput_t>(comm, 10, 2, 1, false);
92 doTest<idInput_t>(comm, 10, 2, 1, true);
93 doTest<idInput_t>(comm, 10, 1, 2, true);
94 doTest<idInput_t>(comm, 10, 1, 2, false);
95 doTest<idInput_t>(comm, 10, 1, -1, false);
96 doTest<idInput_t>(comm, 10, 1, -1, true);
97 doTest<idInput_t>(comm, 10, 2, -1, false);
98
99 if(bCanTestDegreeAsWeights) {
100 doTest<idInput_t>(comm, 10, 1, 1, true, true); // with degreeAsWeights
101 }
102}
103
104int main(int narg, char *arg[])
105{
106 Tpetra::ScopeGuard tscope(&narg, &arg);
107 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
108
109 int rank = comm->getRank();
110
111 // do some tests with BasicIdentifierAdapter
112 runTestSuite<basic_idInput_t>(comm, false);
113
114 // now some tests with XpetraCrsGraphAdapter
115 // Note that right now these are all going to produce the same graph
116 // metrics but could be developed further
117 runTestSuite<graph_idInput_t>(comm, true);
118
119 comm->barrier();
120 if (rank==0)
121 std::cout << "PASS" << std::endl;
122}
123
124// to validate the results, we call evaluate_adapter_results which is
125// templated so it can, for example, check graph metrics only for the graph
126// adapter. Currently both basic and graph adapter setup imbalance metrics so
127// we also do a check on that with a universal call to
128// evaluate_imbalance_results. Currently this needs no specialization.
129// If we add more adapters later this could be flexible to accomodate that.
130template<class idInput_t>
131void evaluate_imbalance_results(RCP<const Comm<int> > comm,
132 RCP<Zoltan2::EvaluatePartition<idInput_t>> metricObject, int numLocalObj,
133 int nWeights, int original_numLocalParts, bool givePartSizes) {
134 int fail = 0;
135 int rank = comm->getRank();
136
137 zscalar_t object_count_imbalance;
138 try{
139 object_count_imbalance = metricObject->getObjectCountImbalance();
140 if(rank == 0) {
141 cout << "Object imbalance: " << object_count_imbalance << endl;
142 }
143 }
144 catch (std::exception &e){
145 fail=1;
146 }
147 TEST_FAIL_AND_EXIT(*comm, fail==0, "getObjectCountImbalance", 1);
148
149 if (nWeights > 0){
150 try{
151 for (int i=0; i < nWeights; i++){
152 zscalar_t imb = metricObject->getWeightImbalance(i);
153 if(rank == 0){
154 cout << "Weight " << i << " imbalance: " << imb << endl;
155 }
156 }
157 }
158 catch (std::exception &e){
159 fail=10;
160 }
161 if (!fail && nWeights > 1){
162 try{
163 zscalar_t imb = metricObject->getNormedImbalance();
164 if(rank == 0){
165 cout << "Normed weight imbalance: " << imb << endl;
166 }
167 }
168 catch (std::exception &e){
169 fail=11;
170 }
171 }
172 }
173 TEST_FAIL_AND_EXIT(*comm, fail==0, "get imbalances", 1);
174}
175
176template<class idInput_t>
177void evaluate_adapter_results(RCP<const Comm<int> > comm,
178 RCP<Zoltan2::EvaluatePartition<idInput_t>> metricObject, int numLocalObj,
179 int nWeights, int original_numLocalParts, bool givePartSizes) {
180 throw std::logic_error("evaluate_result not implemented.");
181}
182
183template<>
184void evaluate_adapter_results<graph_idInput_t>(RCP<const Comm<int> > comm,
185 RCP<Zoltan2::EvaluatePartition<graph_idInput_t>> metricObject, int numLocalObj,
186 int nWeights, int original_numLocalParts, bool givePartSizes) {
187 int fail = 0;
188 int rank = comm->getRank();
189
190 int total_edge_cut = -1;
191 try{
192 // TODO: the unweighted getTotalEdgeCut is an integer
193 // maybe the API should be changed for this and other similar cases
194 total_edge_cut = static_cast<int>(metricObject->getTotalEdgeCut());
195 if(rank == 0){
196 cout << "Total Edge Cut: " << total_edge_cut << endl;
197 }
198 }
199 catch (std::exception &e){
200 fail=1;
201 }
202 TEST_FAIL_AND_EXIT(*comm, fail==0, "getTotalEdgeCut", 1);
203
204 int max_edge_cut = -1;
205 try{
206 max_edge_cut = static_cast<int>(metricObject->getMaxEdgeCut());
207 if(rank == 0){
208 cout << "Max Edge Cut: " << max_edge_cut << endl;
209 }
210 }
211 catch (std::exception &e){
212 fail=1;
213 }
214 TEST_FAIL_AND_EXIT(*comm, fail==0, "getMaxEdgeCut", 1);
215
216 int total_messages = -1;
217 try{
218 total_messages = static_cast<int>(metricObject->getTotalMessages());
219 if(rank == 0){
220 cout << "Total Messages: " << total_messages << endl;
221 }
222 }
223 catch (std::exception &e){
224 fail=1;
225 }
226 TEST_FAIL_AND_EXIT(*comm, fail==0, "getTotalMessages", 1);
227
228 int max_messages = -1;
229 try{
230 max_messages = static_cast<int>(metricObject->getMaxMessages());
231 if(rank == 0){
232 cout << "Max Messages: " << max_messages << endl;
233 }
234 }
235 catch (std::exception &e){
236 fail=1;
237 }
238 TEST_FAIL_AND_EXIT(*comm, fail==0, "getMaxMessages", 1);
239
240 // Now let's check our numbers.
241 // Here we do a calculation of what getTotalEdgeCut should return based on
242 // how we set things up in create_adapter.
243 // Currently the algorithm simply has every object create two links, one
244 // to the first global id and one to the last.
245 // Two of the procs will contain one of those global ids so they only have
246 // edge cuts equal to numLocalObjs to send to the other.
247 // So that is the (2 * numLocalObjs) term.
248 // All other procs will contain neither of those global ids so they have
249 // to send their objects to two procs.
250 // So that is the ((num_procs-2) * numLocalObjs * 2 term.
251 int num_procs = comm->getSize();
252 int expected_total_edge_cuts = (num_procs == 1) ? 0 :
253 (2 * numLocalObj) + ((num_procs-2) * numLocalObj * 2);
254 TEST_FAIL_AND_EXIT(*comm, total_edge_cut == expected_total_edge_cuts,
255 "getTotalEdgeCut is not the expected ", 1);
256
257 // we can also calculate max edge cuts
258 // if num_procs 1, then it's 0
259 // if num_procs 2, then it's numLocalObjs
260 // otherwise it's 2 * numLocalObjs because at least one proc is sending
261 // to two other procs
262 int expected_max_edge_cuts = (num_procs == 1) ? 0 :
263 (num_procs == 2) ? numLocalObj : numLocalObj * 2;
264 TEST_FAIL_AND_EXIT(*comm, max_edge_cut == expected_max_edge_cuts,
265 "getMaxEdgeCut is not the expected value", 1);
266
267 // now check total messages - in present form we can simply divide but in
268 // future things could be generalized
269 int expected_total_messages = expected_total_edge_cuts / numLocalObj;
270 TEST_FAIL_AND_EXIT(*comm, total_messages == expected_total_messages,
271 "getTotalMessages is not the expected value", 1);
272
273 // now check max messages - in present form we can simply divide but in
274 // future things could be more generalized
275 int expected_max_messages = expected_max_edge_cuts / numLocalObj;
276 TEST_FAIL_AND_EXIT(*comm, max_messages == expected_max_messages,
277 "getMaxMessages is not the expected value", 1);
278
279 evaluate_imbalance_results(comm, metricObject,
280 numLocalObj, nWeights, original_numLocalParts, givePartSizes);
281}
282
283// for basic_idInput_t we just call the common evaluate_imbalance_results
284// no other specialized data to consider
285template<>
286void evaluate_adapter_results<basic_idInput_t>(RCP<const Comm<int> > comm,
287 RCP<Zoltan2::EvaluatePartition<basic_idInput_t>> metricObject, int numLocalObj,
288 int nWeights, int original_numLocalParts, bool givePartSizes) {
289 evaluate_imbalance_results(comm, metricObject,
290 numLocalObj, nWeights, original_numLocalParts, givePartSizes);
291}
292
293template<class idInput_t>
294idInput_t * create_adapter(RCP<const Comm<int> > comm,
295 int numLocalObj, zgno_t *myGids,
296 std::vector<const zscalar_t *> & weights,
297 std::vector<int> & strides,
298 bool useDegreeAsWeight) {
299 throw std::logic_error("create_adapter not implemented.");
300}
301
302template<>
304 int numLocalObj, zgno_t *myGids,
305 std::vector<const zscalar_t *> & weights,
306 std::vector<int> & strides,
307 bool useDegreeAsWeight) {
308
309 typedef Tpetra::Map<zlno_t, zgno_t> map_t;
310 typedef Tpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t> matrix_t;
311
312 const zgno_t gNvtx = numLocalObj * comm->getSize();
313 const Teuchos::ArrayView<const zgno_t> indexList(myGids, numLocalObj);
314 Teuchos::RCP<const map_t> map = rcp(new map_t(gNvtx, indexList, 0, comm));
315
316 // Make some stuff in the graph
317 size_t maxRowLen = 2;
318 Teuchos::RCP<matrix_t> matrix = rcp(new matrix_t(map, maxRowLen));
319
320 // I picked this graph as a simple test case.
321 // Something we can easily calculate the final result for as we'd like to
322 // validate this but not end up rewriting the algorithm we are testing.
323 // I have each graph element create two links to the
324 // first global index and last. That means two procs will have edge cuts
325 // equal to their numLocalObj while the rest will have 2 * numLocalObj
326 //
327 // Two of the procs will have only 1 message to send
328 // The other procs will have 2 messages to send
329 // Message max is 2
330 // Message total is going to be (2)*2 + (numProcs-2)*2
331 Teuchos::Array<zgno_t> col(2);
332 Teuchos::Array<zscalar_t> val(2); val[0] = 1.; val[1] = 1.;
333 zgno_t first_id = map->getMinAllGlobalIndex();
334 zgno_t last_id = map->getMaxAllGlobalIndex();
335 for (zlno_t i = 0; i < numLocalObj; i++) {
336 zgno_t id = map->getGlobalElement(i);
337 col[0] = first_id;
338 col[1] = last_id;
339 matrix->insertGlobalValues(id, col(), val());
340 }
341
342 matrix->fillComplete(map, map);
343
344 size_t nVwgts = weights.size();
345 graph_idInput_t * adapter = new graph_idInput_t(matrix->getCrsGraph(), nVwgts);
346
347 // Set the weights
348 for (size_t j = 0; j < nVwgts; j++) {
349 adapter->setWeights(weights[j], 1, j);
350 }
351
352 // set degreeAsWeight if enabled
353 if(useDegreeAsWeight) {
354 for (size_t j = 0; j < nVwgts; j++) {
355 adapter->setWeightIsDegree(j);
356 }
357 }
358
359 return adapter;
360}
361
362template<>
364 int numLocalObj, zgno_t *myGids,
365 std::vector<const zscalar_t *> & weights,
366 std::vector<int> & strides,
367 bool useDegreeAsWeight) {
368 // useDegreeAsWeight is ignored
369 return new basic_idInput_t(numLocalObj, myGids, weights, strides);
370}
371
372// Assumes numLocalObj is the same on every process.
373template<class idInput_t>
374void doTest(RCP<const Comm<int> > comm, int numLocalObj,
375 int nWeights, int numLocalParts, bool givePartSizes, bool useDegreeAsWeight)
376{
378
379 typedef typename idInput_t::part_t part_t;
380
381 int rank = comm->getRank();
382
383 int original_numLocalParts = numLocalParts; // save for log and error checking
384
385 int nprocs = comm->getSize();
386 int fail=0;
387 srand(rank+1);
388 bool testEmptyParts = (numLocalParts < 1);
389 int numGlobalParts = 0;
390
391 if (testEmptyParts){
392 numGlobalParts = nprocs / 2;
393 if (numGlobalParts >= 1)
394 numLocalParts = (rank < numGlobalParts ? 1 : 0);
395 else{
396 numLocalParts = 1;
397 testEmptyParts = false;
398 }
399 }
400 else{
401 numGlobalParts = nprocs * numLocalParts;
402 }
403
404 if(rank == 0) {
405 cout << endl
406 << "Test: number of weights " << nWeights
407 << ", desired number of parts " << numGlobalParts
408 << ", calculated num local parts " << numLocalParts
409 << ", original num local parts " << original_numLocalParts
410 << (givePartSizes ? ", with differing part sizes." :
411 ", with uniform part sizes.")
412 << ", Number of procs " << nprocs
413 << ", each with " << numLocalObj << " objects, part = rank."
414 << (useDegreeAsWeight ? ", use degree as weights" : "")
415 << endl;
416 }
417
418 // An environment. This is usually created by the problem.
419
420 Teuchos::ParameterList pl("test list");
421 pl.set("num_local_parts", numLocalParts);
422
423 RCP<const Zoltan2::Environment> env =
424 rcp(new Zoltan2::Environment(pl, comm));
425
426 // A simple identifier map. Usually created by the model.
427
428 zgno_t *myGids = new zgno_t [numLocalObj];
429 for (int i=0, x=rank*numLocalObj; i < numLocalObj; i++, x++){
430 myGids[i] = x;
431 }
432
433 // Part sizes. Usually supplied by the user to the Problem.
434 // Then the problem supplies them to the Solution.
435
436 int partSizeDim = (givePartSizes ? (nWeights ? nWeights : 1) : 0);
437 ArrayRCP<ArrayRCP<part_t> > ids(partSizeDim);
438 ArrayRCP<ArrayRCP<zscalar_t> > sizes(partSizeDim);
439
440 if (givePartSizes && numLocalParts > 0){
441 part_t *myParts = new part_t [numLocalParts];
442 myParts[0] = rank * numLocalParts;
443 for (int i=1; i < numLocalParts; i++)
444 myParts[i] = myParts[i-1] + 1;
445 ArrayRCP<part_t> partNums(myParts, 0, numLocalParts, true);
446
447 zscalar_t sizeFactor = nprocs/2 - rank;
448 if (sizeFactor < 0) sizeFactor *= -1;
449 sizeFactor += 1;
450
451 for (int dim=0; dim < partSizeDim; dim++){
452 zscalar_t *psizes = new zscalar_t [numLocalParts];
453 for (int i=0; i < numLocalParts; i++)
454 psizes[i] = sizeFactor;
455 sizes[dim] = arcp(psizes, 0, numLocalParts, true);
456 ids[dim] = partNums;
457 }
458 }
459
460 // An input adapter with random weights. Created by the user.
461
462 std::vector<const zscalar_t *> weights;
463 std::vector<int> strides; // default to 1
464
465 int len = numLocalObj*nWeights;
466 ArrayRCP<zscalar_t> wgtBuf;
467 zscalar_t *wgts = NULL;
468
469 if (len > 0){
470 wgts = new zscalar_t [len];
471 wgtBuf = arcp(wgts, 0, len, true);
472 for (int i=0; i < len; i++)
473 wgts[i] = (zscalar_t(rand()) / zscalar_t(RAND_MAX)) + 1.0;
474 }
475
476 for (int i=0; i < nWeights; i++, wgts+=numLocalObj)
477 weights.push_back(wgts);
478
479 idInput_t *ia = NULL;
480
481 try {
482 ia = create_adapter<idInput_t>(comm, numLocalObj, myGids, weights, strides, useDegreeAsWeight);
483 }
484 catch (std::exception &e){
485 fail=1;
486 }
487
488 TEST_FAIL_AND_EXIT(*comm, fail==0, "create adapter", 1);
489
490 // A solution (usually created by a problem)
491
492 RCP<Zoltan2::PartitioningSolution<idInput_t> > solution;
493
494 try{
495 if (givePartSizes)
497 env, comm, nWeights,
498 ids.view(0,partSizeDim), sizes.view(0,partSizeDim)));
499 else
501 env, comm, nWeights));
502 }
503 catch (std::exception &e){
504 fail=1;
505 }
506
507 TEST_FAIL_AND_EXIT(*comm, fail==0, "create solution", 1);
508
509 // Part assignment for my objects: The algorithm usually calls this.
510
511 part_t *partNum = new part_t [numLocalObj];
512 ArrayRCP<part_t> partAssignment(partNum, 0, numLocalObj, true);
513 for (int i=0; i < numLocalObj; i++)
514 partNum[i] = rank;
515
516 solution->setParts(partAssignment);
517
518 // create metric object (also usually created by a problem)
519
520 RCP<quality_t> metricObject;
521
522 try{
523 metricObject = rcp(new quality_t(ia, &pl, comm, solution.getRawPtr()));
524 }
525 catch (std::exception &e){
526 fail=1;
527 }
528 TEST_FAIL_AND_EXIT(*comm, fail==0, "compute metrics", 1);
529
530 try{
531 if(rank == 0){
532 metricObject->printMetrics(cout);
533 }
534 }
535 catch (std::exception &e){
536 fail=1;
537 }
538 TEST_FAIL_AND_EXIT(*comm, fail==0, "print metrics", 1);
539
540 // will call TEST_FAIL_AND_EXIT at each internal step
541 evaluate_adapter_results<idInput_t>(comm, metricObject,
542 numLocalObj, nWeights, original_numLocalParts, givePartSizes);
543
544 delete ia;
545}
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
Zoltan2::BasicUserTypes< zscalar_t, zlno_t, zgno_t > user_t
Definition: Metric.cpp:74
Zoltan2::BasicIdentifierAdapter< user_t > basic_idInput_t
Definition: Metric.cpp:77
Zoltan2::XpetraCrsGraphAdapter< tcrsGraph_t, user_t > graph_idInput_t
Definition: Metric.cpp:81
idInput_t * create_adapter(RCP< const Comm< int > > comm, int numLocalObj, zgno_t *myGids, std::vector< const zscalar_t * > &weights, std::vector< int > &strides, bool useDegreeAsWeight)
Definition: Metric.cpp:294
void evaluate_adapter_results< graph_idInput_t >(RCP< const Comm< int > > comm, RCP< Zoltan2::EvaluatePartition< graph_idInput_t > > metricObject, int numLocalObj, int nWeights, int original_numLocalParts, bool givePartSizes)
Definition: Metric.cpp:184
void doTest(RCP< const Comm< int > > comm, int numLocalObj, int nWeights, int numLocalParts, bool givePartSizes, bool useDegreeAsWeight=false)
Definition: Metric.cpp:374
basic_idInput_t * create_adapter< basic_idInput_t >(RCP< const Comm< int > > comm, int numLocalObj, zgno_t *myGids, std::vector< const zscalar_t * > &weights, std::vector< int > &strides, bool useDegreeAsWeight)
Definition: Metric.cpp:363
void evaluate_imbalance_results(RCP< const Comm< int > > comm, RCP< Zoltan2::EvaluatePartition< idInput_t > > metricObject, int numLocalObj, int nWeights, int original_numLocalParts, bool givePartSizes)
Definition: Metric.cpp:131
graph_idInput_t * create_adapter< graph_idInput_t >(RCP< const Comm< int > > comm, int numLocalObj, zgno_t *myGids, std::vector< const zscalar_t * > &weights, std::vector< int > &strides, bool useDegreeAsWeight)
Definition: Metric.cpp:303
void evaluate_adapter_results(RCP< const Comm< int > > comm, RCP< Zoltan2::EvaluatePartition< idInput_t > > metricObject, int numLocalObj, int nWeights, int original_numLocalParts, bool givePartSizes)
Definition: Metric.cpp:177
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > tcrsGraph_t
Definition: Metric.cpp:80
void evaluate_adapter_results< basic_idInput_t >(RCP< const Comm< int > > comm, RCP< Zoltan2::EvaluatePartition< basic_idInput_t > > metricObject, int numLocalObj, int nWeights, int original_numLocalParts, bool givePartSizes)
Definition: Metric.cpp:286
void runTestSuite(RCP< const Comm< int > > comm, bool bCanTestDegreeAsWeights)
Definition: Metric.cpp:85
Defines the BasicIdentifierAdapter class.
Defines the EvaluatePartition class.
common code used by tests
float zscalar_t
Tpetra::Map ::local_ordinal_type zlno_t
Tpetra::Map ::global_ordinal_type zgno_t
Defines XpetraCrsGraphAdapter class.
int main()
This class represents a collection of global Identifiers and their associated weights,...
A simple class that can be the User template argument for an InputAdapter.
The user parameters, debug, timing and memory profiling output objects, and error checking methods.
A class that computes and returns quality metrics.
A PartitioningSolution is a solution to a partitioning problem.
Provides access for Zoltan2 to Xpetra::CrsGraph data.
void setWeightIsDegree(int idx)
Specify an index for which the weight should be the degree of the entity.
void setWeights(const scalar_t *val, int stride, int idx)
Provide a pointer to weights for the primary entity type.
static const std::string fail
Tpetra::Map map_t
Definition: mapRemotes.cpp:16
static ArrayRCP< ArrayRCP< zscalar_t > > weights
SparseMatrixAdapter_t::part_t part_t
Zoltan2::EvaluatePartition< matrixAdapter_t > quality_t