Zoltan2
Loading...
Searching...
No Matches
coloring1.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#include <iostream>
46#include <fstream>
47#include <limits>
48#include <vector>
49#include <Teuchos_RCP.hpp>
50#include <Teuchos_ParameterList.hpp>
51#include <Teuchos_CommandLineProcessor.hpp>
52#include <Tpetra_CrsMatrix.hpp>
53#include <Tpetra_Vector.hpp>
54#include <MatrixMarket_Tpetra.hpp>
58
59using Teuchos::RCP;
60
62// Program to demonstrate use of Zoltan2 to color a TPetra matrix
63// (read from a MatrixMarket file or generated by Galeri::Xpetra).
64// We assume the matrix is structurally symmetric.
65// Usage:
66// a.out [--inputFile=filename] [--outputFile=outfile] [--verbose]
67// [--x=#] [--y=#] [--z=#] [--matrix={Laplace1D,Laplace2D,Laplace3D}
69
71// Eventually want to use Teuchos unit tests to vary z2TestLO and
72// GO. For now, we set them at compile time based on whether Tpetra
73// is built with explicit instantiation on. (in Zoltan2_TestHelpers.hpp)
74
78
79typedef Tpetra::CrsMatrix<z2TestScalar, z2TestLO, z2TestGO> SparseMatrix;
80typedef Tpetra::Vector<z2TestScalar, z2TestLO, z2TestGO> Vector;
81typedef Vector::node_type Node;
82
83typedef Tpetra::Import<z2TestLO, z2TestGO> Import;
84
86
87int validateColoring(RCP<SparseMatrix> A, int *color)
88// returns 0 if coloring is valid, nonzero if invalid
89{
90 if (A->getRowMap()->getComm()->getRank() == 0)
91 std::cout << "validate coloring: local graph, distance-one" << std::endl;
92 int nconflicts = 0;
93 typename tcrsMatrix_t::local_inds_host_view_type indices;
94 typename tcrsMatrix_t::values_host_view_type values;
95 // Count conflicts in the graph.
96 // Loop over local rows, treat local column indices as edges.
97 zlno_t n = A->getLocalNumRows();
98 for (zlno_t i=0; i<n; i++) {
99 A->getLocalRowView(i, indices, values);
100 for (zlno_t j=0; j<static_cast<zlno_t>(indices.extent(0)); j++) {
101 if ((indices[j]<n) && (indices[j]!=i) && (color[i]==color[indices[j]])){
102 nconflicts++;
103 //std::cout << "Debug: found conflict (" << i << ", " << indices[j] << ")" << std::endl;
104 }
105 }
106 }
107
108 return nconflicts;
109}
110
111int validateDistributedColoring(RCP<SparseMatrix> A, int *color){
112 if (A->getRowMap()->getComm()->getRank() == 0)
113 std::cout << "validate coloring: distributed, distance-one" << std::endl;
114 int nconflicts = 0;
115 RCP<const SparseMatrix::map_type> rowMap = A->getRowMap();
116 RCP<const SparseMatrix::map_type> colMap = A->getColMap();
117 Vector R = Vector(rowMap);
118 //put colors in the scalar entries of R
119 for(size_t i = 0; i < A->getLocalNumRows(); i++){
120 R.replaceLocalValue(i, color[i]);
121 }
122
123 Vector C = Vector(colMap);
124 Import imp = Import(rowMap, colMap);
125 C.doImport(R, imp, Tpetra::REPLACE);
126
127 typename tcrsMatrix_t::local_inds_host_view_type indices;
128 typename tcrsMatrix_t::values_host_view_type values;
129
130 //count conflicts in the graph
131 //loop over local rows, treat local column indices as edges
132 size_t n = A->getLocalNumRows();
133 auto colorData = C.getData();
134 for(size_t i = 0; i < n; i++){
135 A->getLocalRowView(i, indices, values);
136 for(size_t j = 0; j < indices.extent(0); j++){
137 if(values[j] == 0) continue; //this catches removed entries.
138 if( (rowMap->getGlobalElement(i) != colMap->getGlobalElement(indices[j])) && (color[i] == colorData[indices[j]]) ){
139 nconflicts++;
140 }
141 }
142 }
143
144 return nconflicts;
145}
146
147int validateDistributedDistance2Coloring(RCP<SparseMatrix> A, int* color){
148
149 if (A->getRowMap()->getComm()->getRank() == 0)
150 std::cout << "validate coloring: distributed, distance-two" << std::endl;
151
152 //To check distance-2 conflicts, we square the input matrix and check
153 //for distance-1 conflicts on the squared matrix.
154 RCP<SparseMatrix> S = rcp(new SparseMatrix(A->getRowMap(), 0));
155 Tpetra::MatrixMatrix::Multiply(*A, true, *A, false, *S);
156
157 return validateDistributedColoring(S,color);
158}
159
160int checkBalance(zlno_t n, int *color)
161// Check size of color classes
162{
163 // Find max color
164 int maxColor = 0;
165 for (zlno_t i=0; i<n; i++) {
166 if (color[i] > maxColor) maxColor = color[i];
167 }
168
169 // Compute color class sizes
170 Teuchos::Array<int> colorCount(maxColor+1);
171 for (zlno_t i=0; i<n; i++) {
172 colorCount[color[i]]++;
173 }
174
175 // Find min and max, excluding color 0.
176 int smallest = 1;
177 int largest = 1;
178 zlno_t small = colorCount[1];
179 zlno_t large = colorCount[1];
180 for (int i=1; i<=maxColor; i++){
181 if (colorCount[i] < small){
182 small = colorCount[i];
183 smallest = i;
184 }
185 if (colorCount[i] > large){
186 large = colorCount[i];
187 largest = i;
188 }
189 }
190
191 //std::cout << "Color size[0:2] = " << colorCount[0] << ", " << colorCount[1] << ", " << colorCount[2] << std::endl;
192 std::cout << "Largest color class = " << largest << " with " << colorCount[largest] << " vertices." << std::endl;
193 std::cout << "Smallest color class = " << smallest << " with " << colorCount[smallest] << " vertices." << std::endl;
194
195 return 0;
196}
197
199int main(int narg, char** arg)
200{
201 std::string inputFile = ""; // Matrix Market file to read
202 std::string outputFile = ""; // Output file to write
203 std::string colorAlg = "SerialGreedy"; // Default algorithm is the serial greedy
204 bool verbose = false; // Verbosity of output
205 bool timing = false; // If true, report coloring times.
206 int testReturn = 0;
207 bool recolorDegrees = false;
208 std::string prepartition = ""; // Call Zoltan2 partitioning to better distribute
209 // the graph before coloring
210 bool prepartition_rows = false; // True if we're prepartitioning wrt rows
211 bool prepartition_nonzeros = false; // True if prepartitioning wrt matrix nonzeros
212
213 int localColors = 0;
214 int totalColors = 0;
216 Tpetra::ScopeGuard tscope(&narg, &arg);
217 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
218 int me = comm->getRank();
219 int serialThreshold = 0;
220 // Read run-time options.
221 Teuchos::CommandLineProcessor cmdp (false, false);
222 cmdp.setOption("colorMethod", &colorAlg,
223 "Coloring algorithms supported: SerialGreedy, D1, D1-2GL, D2, PD2");
224 cmdp.setOption("inputFile", &inputFile,
225 "Name of a Matrix Market file in the data directory; "
226 "if not specified, a matrix will be generated by Galeri.");
227 cmdp.setOption("outputFile", &outputFile,
228 "Name of file to write the coloring");
229 cmdp.setOption("verbose", "quiet", &verbose,
230 "Print messages and results.");
231 cmdp.setOption("prepartition", &prepartition,
232 "Partition the input graph for better initial distribution;"
233 "valid values are rows and nonzeros");
234 cmdp.setOption("serialThreshold", &serialThreshold,
235 "number of vertices to recolor in serial");
236 cmdp.setOption("recolorDegrees","recolorRandom",&recolorDegrees,
237 "recolor based on vertex degrees or random numbers");
238 cmdp.setOption("timing", "notimes", &timing,
239 "report how long coloring takes");
240 std::cout << "Starting everything" << std::endl;
241
243 // Even with cmdp option "true", I get errors for having these
244 // arguments on the command line. (On redsky build)
245 // KDDKDD Should just be warnings, right? Code should still work with these
246 // KDDKDD params in the create-a-matrix file. Better to have them where
247 // KDDKDD they are used.
248 int xdim=10;
249 int ydim=10;
250 int zdim=10;
251 std::string matrixType("Laplace3D");
252
253 cmdp.setOption("x", &xdim,
254 "number of gridpoints in X dimension for "
255 "mesh used to generate matrix.");
256 cmdp.setOption("y", &ydim,
257 "number of gridpoints in Y dimension for "
258 "mesh used to generate matrix.");
259 cmdp.setOption("z", &zdim,
260 "number of gridpoints in Z dimension for "
261 "mesh used to generate matrix.");
262 cmdp.setOption("matrix", &matrixType,
263 "Matrix type: Laplace1D, Laplace2D, or Laplace3D");
264
266 // Coloring options to test.
268 std::string colorMethod("FirstFit");
269 //int balanceColors = 0;
270 cmdp.setOption("color_choice", &colorMethod,
271 "Color choice method: FirstFit, LeastUsed, Random, RandomFast");
272 // cmdp.setOption("balance_colors", &balanceColors,
273 // "Balance the size of color classes: 0/1 for false/true");
274
276 cmdp.parse(narg, arg);
277
278 if(prepartition != ""){
279 if(prepartition == "rows") prepartition_rows = true;
280 else if (prepartition == "nonzeros") prepartition_nonzeros = true;
281 else {
282 std::cout << "Invalid value of prepartition option " << prepartition
283 << std::endl;
284 std::cout << "No prepartitioning will be done" <<std::endl;
285 }
286 }
287
288 RCP<UserInputForTests> uinput;
289
290 if (inputFile != ""){ // Input file specified; read a matrix
291 uinput = rcp(new UserInputForTests(testDataFilePath, inputFile, comm, true));
292 }
293 else // Let Galeri generate a matrix
294
295 uinput = rcp(new UserInputForTests(xdim, ydim, zdim, matrixType, comm, true, true));
296
297 RCP<SparseMatrix> Matrix = uinput->getUITpetraCrsMatrix();
298
299 if (me == 0)
300 std::cout << "NumRows = " << Matrix->getGlobalNumRows() << std::endl
301 << "NumNonzeros = " << Matrix->getGlobalNumEntries() << std::endl
302 << "NumProcs = " << comm->getSize() << std::endl;
303 if(prepartition_rows || prepartition_nonzeros){
304 std::cout<<comm->getRank() <<": Starting to pre-partition, creating adapter\n";
305 //compute new partition of matrix
306 std::unique_ptr<SparseMatrixAdapter> zadapter;
307 if(prepartition_nonzeros){
308 zadapter = std::unique_ptr<SparseMatrixAdapter>(new SparseMatrixAdapter(Matrix, 1));
309 zadapter->setRowWeightIsNumberOfNonZeros(0);
310 } else {
311 zadapter = std::unique_ptr<SparseMatrixAdapter>(new SparseMatrixAdapter(Matrix));
312 }
313
314 std::cout<<comm->getRank()<<": created adapter, creating PartitioningProblem\n";
315 Teuchos::ParameterList zparams;
316 zparams.set("algorithm","parmetis");
317 zparams.set("imbalance_tolerance", 1.05);
318 zparams.set("partitioning_approach","partition");
320 zproblem(zadapter.get(), &zparams);
321 std::cout<<comm->getRank()<<": created PartitioningProblem, starting to solve\n";
322 zproblem.solve();
323 std::cout<<comm->getRank()<<": solved Partitioning Problem\n";
324 //print partition characteristics before and after
325 std::cout<<comm->getRank()<<": applying partition\n";
327 quality_t evalbef(zadapter.get(), &zparams, comm, NULL);
328 if(me == 0){
329 std::cout<<"BEFORE PREPARTITION: Partition statistics:" << std::endl;
330 evalbef.printMetrics(std::cout);
331 }
332
333 quality_t evalaft(zadapter.get(), &zparams, comm, &zproblem.getSolution());
334 if(me == 0){
335 std::cout<<"AFTER PREPARTITION: Partition statistics:"<<std::endl;
336 evalaft.printMetrics(std::cout);
337 }
338 std::cout<<comm->getRank()<<": done evaluating, migrating matrix to use new partitioning\n";
339 //Migrate matrix to the new partition
340 RCP<SparseMatrix> newMatrix;
341 zadapter->applyPartitioningSolution(*Matrix, newMatrix, zproblem.getSolution());
342
343 std::cout<<comm->getRank()<<": done applying, replacing old matrix with new one\n";
344 Matrix = newMatrix;
345 std::cout<<comm->getRank()<<": done replacing, finished partitioning\n";
346 }
348 Teuchos::ParameterList params;
349 params.set("color_choice", colorMethod);
350 params.set("color_method", colorAlg);
351 params.set("verbose", verbose);
352 params.set("timing", timing);
353 params.set("serial_threshold",serialThreshold);
354 params.set("recolor_degrees",recolorDegrees);
355
356 //params.set("balance_colors", balanceColors); // TODO
357
359 SparseMatrixAdapter adapter(Matrix);
360
362 try
363 {
364 Zoltan2::ColoringProblem<SparseMatrixAdapter> problem(&adapter, &params);
365 std::cout << "Going to color " << colorAlg << std::endl;
366 problem.solve();
367
369 size_t checkLength;
370 int *checkColoring = nullptr;
371
373
374 std::cout << "Going to get results" << std::endl;
375 // Check that the solution is really a coloring
376 checkLength = soln->getColorsSize();
377 if(checkLength > 0){
378 checkColoring = soln->getColors();
379 localColors = soln->getNumColors();
380 }
381
382 //count the number of colors used globally
383 Teuchos::reduceAll<int,int>(*comm, Teuchos::REDUCE_MAX, 1, &localColors, &totalColors);
384 if (outputFile != "") {
385 std::ofstream colorFile;
386
387 // Write coloring to file,
388 // each process writes local coloring to a separate file
389 //std::string fname = outputFile + "." + me;
390 std::stringstream fname;
391 fname << outputFile << "." << comm->getSize() << "." << me;
392 colorFile.open(fname.str().c_str());
393 for (size_t i=0; i<checkLength; i++){
394 colorFile << " " << checkColoring[i] << std::endl;
395 }
396 colorFile.close();
397 }
398
399 // Print # of colors on each proc.
400 std::cout << "No. of colors on proc " << me << " : " << localColors << std::endl;
401
402 std::cout << "Going to validate the soln" << std::endl;
403 auto rowInds = Matrix->getRowMap()->getMyGlobalIndices();
404 Matrix->resumeFill();
405
406 //We use a squared matrix to validated distance-2 colorings.
407 //If any diagonals are present in the input matrix before we
408 //square the matrix, the resultant matrix may have both distance-1
409 //and distance-2 paths. By removing the diagonals before squaring,
410 //we ensure that our distance-2 validation only checks distance-2
411 //paths. This makes it suitable for validating partial distance-2 colorings.
412
413 std::cout<<"Got row indices, replacing diagonals\n";
414 //loop through all rows
415 for(size_t i = 0; i < rowInds.size(); i++){
416 typename Kokkos::View<zgno_t*>::HostMirror idx("idx",1);
417 typename Kokkos::View<zscalar_t*>::HostMirror val("val",1);
418 Kokkos::deep_copy(idx,rowInds[i]);
419 Kokkos::deep_copy(val, 0.);
420 //get the entries in the current row
421 size_t numEntries = Matrix->getNumEntriesInGlobalRow(rowInds[i]);
422 typename tcrsMatrix_t::nonconst_global_inds_host_view_type inds("Indices", numEntries);
423 typename tcrsMatrix_t::nonconst_values_host_view_type vals("Values", numEntries);
424 Matrix->getGlobalRowCopy(rowInds[i], inds, vals, numEntries);
425 //check to see if this row has a diagonal
426 bool hasDiagonal = false;
427 for(size_t j = 0; j < inds.extent(0); j++){
428 if(inds[j] == rowInds[i]) hasDiagonal = true;
429 }
430 //replace the diagonal if it exists.
431 if(hasDiagonal) {
432 if(Matrix->replaceGlobalValues(rowInds[i],idx,val) == Teuchos::OrdinalTraits<zlno_t>::invalid()){
433 std::cout<<"Either !isFillActive() or inds.extent != vals.extent()\n";
434 }
435 //else {
436 // std::cout<<"*******DIAGONAL REPLACED*********\n";
437 //}
438 }
439 }
440 Matrix->fillComplete();
441 // Verify that checkColoring is a coloring
442 if(colorAlg == "D2"){
443 testReturn = validateDistributedColoring(Matrix,checkColoring);
444 testReturn += validateDistributedDistance2Coloring(Matrix, checkColoring);
445 } else if (colorAlg == "PD2"){
446 testReturn = validateDistributedDistance2Coloring(Matrix, checkColoring);
447 }else if(colorAlg == "D1-2GL" || colorAlg == "D1"){
448 testReturn = validateDistributedColoring(Matrix, checkColoring);
449 } else if (checkLength > 0){
450 testReturn = validateColoring(Matrix, checkColoring);
451 }
452
453 // Check balance (not part of pass/fail for now)
454 if(checkLength > 0) checkBalance((zlno_t)checkLength, checkColoring);
455
456 } catch (std::exception &e){
457 std::cout << "Exception caught in coloring" << std::endl;
458 std::cout << e.what() << std::endl;
459 std::cout << "FAIL" << std::endl;
460 return 0;
461 }
462
463 int numGlobalConflicts = 0;
464 Teuchos::reduceAll<int, int>(*comm, Teuchos::REDUCE_MAX, 1, &testReturn, &numGlobalConflicts);
465
466 if (me == 0) {
467 if (numGlobalConflicts > 0){
468 std::cout <<"Number of conflicts found = "<<numGlobalConflicts<<std::endl;
469 std::cout << "Solution is not a valid coloring; FAIL" << std::endl;
470 }else{
471 std::cout << "Used " <<totalColors<<" colors\n";
472 std::cout << "PASS" << std::endl;
473 }
474 }
475
476}
477
Defines the ColoringProblem class.
common code used by tests
float zscalar_t
Tpetra::Map ::local_ordinal_type zlno_t
std::string testDataFilePath(".")
Tpetra::Map ::global_ordinal_type zgno_t
Defines the XpetraCrsMatrixAdapter class.
int main()
ColoringProblem sets up coloring problems for the user.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
ColoringSolution< Adapter > * getSolution()
Get the solution to the problem.
The class containing coloring solution.
int * getColors()
Get (local) color array by raw pointer (no RCP).
size_t getColorsSize()
Get (local) size of color array.
int getNumColors()
Get local number of colors. This is computed from the coloring each time, as this is cheap.
A class that computes and returns quality metrics.
void printMetrics(std::ostream &os) const
Print all metrics.
PartitioningProblem sets up partitioning problems for the user.
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
Provides access for Zoltan2 to Xpetra::CrsMatrix data.
Vector::node_type Node
Definition: coloring1.cpp:81
zlno_t z2TestLO
Definition: coloring1.cpp:75
Tpetra::Import< z2TestLO, z2TestGO > Import
Definition: coloring1.cpp:83
int validateDistributedDistance2Coloring(RCP< SparseMatrix > A, int *color)
Definition: coloring1.cpp:147
int validateDistributedColoring(RCP< SparseMatrix > A, int *color)
Definition: coloring1.cpp:111
Tpetra::CrsMatrix< z2TestScalar, z2TestLO, z2TestGO > SparseMatrix
Definition: coloring1.cpp:79
Zoltan2::XpetraCrsMatrixAdapter< SparseMatrix > SparseMatrixAdapter
Definition: coloring1.cpp:85
Tpetra::Vector< z2TestScalar, z2TestLO, z2TestGO > Vector
Definition: coloring1.cpp:80
zgno_t z2TestGO
Definition: coloring1.cpp:76
int validateColoring(RCP< SparseMatrix > A, int *color)
Definition: coloring1.cpp:87
zscalar_t z2TestScalar
Definition: coloring1.cpp:77
int checkBalance(zlno_t n, int *color)
Definition: coloring1.cpp:160
Zoltan2::EvaluatePartition< matrixAdapter_t > quality_t