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>
79 typedef Tpetra::CrsMatrix<z2TestScalar, z2TestLO, z2TestGO>
SparseMatrix;
80 typedef Tpetra::Vector<z2TestScalar, z2TestLO, z2TestGO>
Vector;
81 typedef Vector::node_type
Node;
83 typedef Tpetra::Import<z2TestLO, z2TestGO>
Import;
90 if (A->getRowMap()->getComm()->getRank() == 0)
91 std::cout <<
"validate coloring: local graph, distance-one" << std::endl;
93 typename tcrsMatrix_t::local_inds_host_view_type indices;
94 typename tcrsMatrix_t::values_host_view_type values;
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]])){
112 if (A->getRowMap()->getComm()->getRank() == 0)
113 std::cout <<
"validate coloring: distributed, distance-one" << std::endl;
115 RCP<const SparseMatrix::map_type> rowMap = A->getRowMap();
116 RCP<const SparseMatrix::map_type> colMap = A->getColMap();
119 for(
size_t i = 0; i < A->getLocalNumRows(); i++){
120 R.replaceLocalValue(i, color[i]);
125 C.doImport(R, imp, Tpetra::REPLACE);
127 typename tcrsMatrix_t::local_inds_host_view_type indices;
128 typename tcrsMatrix_t::values_host_view_type values;
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;
138 if( (rowMap->getGlobalElement(i) != colMap->getGlobalElement(indices[j])) && (color[i] == colorData[indices[j]]) ){
149 if (A->getRowMap()->getComm()->getRank() == 0)
150 std::cout <<
"validate coloring: distributed, distance-two" << std::endl;
154 RCP<SparseMatrix> S = rcp(
new SparseMatrix(A->getRowMap(), 0));
155 Tpetra::MatrixMatrix::Multiply(*A,
true, *A,
false, *S);
165 for (
zlno_t i=0; i<n; i++) {
166 if (color[i] > maxColor) maxColor = color[i];
170 Teuchos::Array<int> colorCount(maxColor+1);
171 for (
zlno_t i=0; i<n; i++) {
172 colorCount[color[i]]++;
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];
185 if (colorCount[i] > large){
186 large = colorCount[i];
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;
201 std::string inputFile =
"";
202 std::string outputFile =
"";
203 std::string colorAlg =
"SerialGreedy";
204 bool verbose =
false;
207 bool recolorDegrees =
false;
208 std::string prepartition =
"";
210 bool prepartition_rows =
false;
211 bool prepartition_nonzeros =
false;
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;
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;
251 std::string matrixType(
"Laplace3D");
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");
268 std::string colorMethod(
"FirstFit");
270 cmdp.setOption(
"color_choice", &colorMethod,
271 "Color choice method: FirstFit, LeastUsed, Random, RandomFast");
276 cmdp.parse(narg, arg);
278 if(prepartition !=
""){
279 if(prepartition ==
"rows") prepartition_rows =
true;
280 else if (prepartition ==
"nonzeros") prepartition_nonzeros =
true;
282 std::cout <<
"Invalid value of prepartition option " << prepartition
284 std::cout <<
"No prepartitioning will be done" <<std::endl;
288 RCP<UserInputForTests> uinput;
290 if (inputFile !=
""){
295 uinput = rcp(
new UserInputForTests(xdim, ydim, zdim, matrixType, comm,
true,
true));
297 RCP<SparseMatrix> Matrix = uinput->getUITpetraCrsMatrix();
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";
306 std::unique_ptr<SparseMatrixAdapter> zadapter;
307 if(prepartition_nonzeros){
309 zadapter->setRowWeightIsNumberOfNonZeros(0);
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";
323 std::cout<<comm->getRank()<<
": solved Partitioning Problem\n";
325 std::cout<<comm->getRank()<<
": applying partition\n";
327 quality_t evalbef(zadapter.get(), &zparams, comm, NULL);
329 std::cout<<
"BEFORE PREPARTITION: Partition statistics:" << std::endl;
330 evalbef.printMetrics(std::cout);
333 quality_t evalaft(zadapter.get(), &zparams, comm, &zproblem.getSolution());
335 std::cout<<
"AFTER PREPARTITION: Partition statistics:"<<std::endl;
336 evalaft.printMetrics(std::cout);
338 std::cout<<comm->getRank()<<
": done evaluating, migrating matrix to use new partitioning\n";
340 RCP<SparseMatrix> newMatrix;
341 zadapter->applyPartitioningSolution(*Matrix, newMatrix, zproblem.getSolution());
343 std::cout<<comm->getRank()<<
": done applying, replacing old matrix with new one\n";
345 std::cout<<comm->getRank()<<
": done replacing, finished partitioning\n";
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);
365 std::cout <<
"Going to color " << colorAlg << std::endl;
370 int *checkColoring =
nullptr;
374 std::cout <<
"Going to get results" << std::endl;
383 Teuchos::reduceAll<int,int>(*comm, Teuchos::REDUCE_MAX, 1, &localColors, &totalColors);
384 if (outputFile !=
"") {
385 std::ofstream colorFile;
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;
400 std::cout <<
"No. of colors on proc " << me <<
" : " << localColors << std::endl;
402 std::cout <<
"Going to validate the soln" << std::endl;
403 auto rowInds = Matrix->getRowMap()->getMyGlobalIndices();
404 Matrix->resumeFill();
413 std::cout<<
"Got row indices, replacing diagonals\n";
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.);
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);
426 bool hasDiagonal =
false;
427 for(
size_t j = 0; j < inds.extent(0); j++){
428 if(inds[j] == rowInds[i]) hasDiagonal =
true;
432 if(Matrix->replaceGlobalValues(rowInds[i],idx,val) == Teuchos::OrdinalTraits<zlno_t>::invalid()){
433 std::cout<<
"Either !isFillActive() or inds.extent != vals.extent()\n";
440 Matrix->fillComplete();
442 if(colorAlg ==
"D2"){
445 }
else if (colorAlg ==
"PD2"){
447 }
else if(colorAlg ==
"D1-2GL" || colorAlg ==
"D1"){
449 }
else if (checkLength > 0){
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;
463 int numGlobalConflicts = 0;
464 Teuchos::reduceAll<int, int>(*comm, Teuchos::REDUCE_MAX, 1, &testReturn, &numGlobalConflicts);
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;
471 std::cout <<
"Used " <<totalColors<<
" colors\n";
472 std::cout <<
"PASS" << std::endl;
int * getColors()
Get (local) color array by raw pointer (no RCP).
ColoringProblem sets up coloring problems for the user.
Defines the ColoringProblem class.
Provides access for Zoltan2 to Xpetra::CrsMatrix data.
Tpetra::CrsMatrix< z2TestScalar, z2TestLO, z2TestGO > SparseMatrix
void solve(bool updateInputData=true)
Direct the problem to create a solution.
Tpetra::Import< z2TestLO, z2TestGO > Import
ColoringSolution< Adapter > * getSolution()
Get the solution to the problem.
int getNumColors()
Get local number of colors. This is computed from the coloring each time, as this is cheap...
Zoltan2::EvaluatePartition< matrixAdapter_t > quality_t
int main(int narg, char **arg)
common code used by tests
Tpetra::Vector< z2TestScalar, z2TestLO, z2TestGO > Vector
Defines the XpetraCrsMatrixAdapter class.
int validateColoring(RCP< SparseMatrix > A, int *color)
int validateDistributedColoring(RCP< SparseMatrix > A, int *color)
Tpetra::Map::local_ordinal_type zlno_t
PartitioningProblem sets up partitioning problems for the user.
int validateDistributedDistance2Coloring(RCP< SparseMatrix > A, int *color)
Zoltan2::XpetraCrsMatrixAdapter< SparseMatrix > SparseMatrixAdapter
A class that computes and returns quality metrics.
size_t getColorsSize()
Get (local) size of color array.
Tpetra::Map::global_ordinal_type zgno_t
void solve(bool updateInputData=true)
Direct the problem to create a solution.
The class containing coloring solution.
int checkBalance(zlno_t n, int *color)
std::string testDataFilePath(".")