52 #include <Teuchos_ParameterList.hpp>
53 #include <Teuchos_RCP.hpp>
54 #include <Teuchos_CommandLineProcessor.hpp>
55 #include <Tpetra_CrsMatrix.hpp>
56 #include <Tpetra_Vector.hpp>
57 #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;
89 typename SparseMatrix::local_inds_host_view_type indices;
90 typename SparseMatrix::values_host_view_type values;
98 for (ii=0; ii<n; ii++) {
99 A->getLocalRowView (ii, indices, values);
100 for (k=0; k< static_cast<z2TestLO>(indices.size()); k++){
104 j = perm[indices[k]];
117 return (bw_left + bw_right + 1);
131 lno_t numRows = origMatrix->getLocalNumRows();
133 std::cout <<
"origMatrix->getLocalNumRows(): " << numRows << std::endl;
136 std::cout <<
"Skipping analysis - matrix is empty" << std::endl;
140 Array<lno_t> oldMatrix(numRows*numRows);
141 Array<lno_t> newMatrix(numRows*numRows);
144 lno_t showSize = numRows;
149 std::cout << std::endl <<
"perm: ";
150 for(lno_t n = 0; n < showSize; ++n) {
151 std::cout <<
" " << perm[n] <<
" ";
153 if(showSize < numRows) {
156 std::cout << std::endl <<
"iperm: ";
157 for(lno_t n = 0; n < showSize; ++n) {
158 std::cout <<
" " << iperm[n] <<
" ";
160 if(showSize < numRows) {
164 std::cout << std::endl;
167 for (lno_t j = 0; j < numRows; ++j) {
168 typename SparseMatrix::local_inds_host_view_type indices;
169 typename SparseMatrix::values_host_view_type wgts;
170 origMatrix->getLocalRowView( j, indices, wgts );
171 for (lno_t n = 0; n < static_cast<lno_t>(indices.size()); ++n) {
172 lno_t i = indices[n];
174 oldMatrix[i + j*numRows] = 1;
175 newMatrix[perm[i] + perm[j]*numRows] = 1;
181 std::cout << std::endl <<
"original graph in matrix form:" << std::endl;
182 for(lno_t y = 0; y < showSize; ++y) {
183 for(lno_t x = 0; x < showSize; ++x) {
184 std::cout <<
" " << oldMatrix[x + y*numRows];
186 if(showSize < numRows) {
189 std::cout << std::endl;
191 if(showSize < numRows) {
192 for(lno_t i = 0; i < showSize; ++i) {
197 std::cout << std::endl;
200 std::cout << std::endl <<
"reordered graph in matrix form:" << std::endl;
201 for(lno_t y = 0; y < showSize; ++y) {
202 for(lno_t x = 0; x < showSize; ++x) {
203 std::cout <<
" " << newMatrix[x + y*numRows];
205 if(showSize < numRows) {
208 std::cout << std::endl;
210 if(showSize < numRows) {
211 for(lno_t i = 0; i < showSize; ++i) {
216 std::cout << std::endl;
221 int mainExecute(
int narg,
char** arg, RCP<
const Teuchos::Comm<int> > comm)
223 std::string inputFile =
"";
224 std::string outputFile =
"";
225 bool verbose =
false;
228 int rank = comm->getRank();
231 Teuchos::CommandLineProcessor cmdp (
false,
false);
232 cmdp.setOption(
"inputFile", &inputFile,
233 "Name of a Matrix Market file in the data directory; "
234 "if not specified, a matrix will be generated by Galeri.");
235 cmdp.setOption(
"outputFile", &outputFile,
236 "Name of file to write the permutation");
237 cmdp.setOption(
"verbose",
"quiet", &verbose,
238 "Print messages and results.");
241 std::cout <<
"Starting everything" << std::endl;
253 std::string matrixType(
"Laplace3D");
255 cmdp.setOption(
"x", &xdim,
256 "number of gridpoints in X dimension for "
257 "mesh used to generate matrix.");
258 cmdp.setOption(
"y", &ydim,
259 "number of gridpoints in Y dimension for "
260 "mesh used to generate matrix.");
261 cmdp.setOption(
"z", &zdim,
262 "number of gridpoints in Z dimension for "
263 "mesh used to generate matrix.");
264 cmdp.setOption(
"matrix", &matrixType,
265 "Matrix type: Laplace1D, Laplace2D, or Laplace3D");
270 std::string orderMethod(
"scotch");
271 cmdp.setOption(
"order_method", &orderMethod,
272 "order_method: natural, random, rcm, scotch");
274 std::string orderMethodType(
"local");
275 cmdp.setOption(
"order_method_type", &orderMethodType,
276 "local or global or both" );
279 cmdp.parse(narg, arg);
282 RCP<UserInputForTests> uinput;
288 if (inputFile !=
""){
296 RCP<SparseMatrix> origMatrix = uinput->getUITpetraCrsMatrix();
299 std::cout <<
"NumRows = " << origMatrix->getGlobalNumRows() << std::endl
300 <<
"NumNonzeros = " << origMatrix->getGlobalNumEntries() << std::endl
301 <<
"NumProcs = " << comm->getSize() << std::endl;
316 Teuchos::ParameterList params;
317 params.set(
"order_method", orderMethod);
318 params.set(
"order_method_type", orderMethodType);
330 std::cout <<
"Going to solve" << std::endl;
341 std::cout <<
"Going to get results" << std::endl;
344 #ifdef MDM // Temp debugging code all of which can be removed for final
345 for(
int checkRank = 0; checkRank < comm->getSize(); ++checkRank ) {
347 if( checkRank == comm->getRank() ) {
348 std::cout <<
"Debug output matrix for rank: " << checkRank << std::endl;
378 if (outputFile !=
"") {
379 std::ofstream permFile;
384 std::stringstream
fname;
385 fname << outputFile <<
"." << comm->getSize() <<
"." << rank;
386 permFile.open(fname.str().c_str());
387 for (
size_t i=0; i<checkLength; i++){
388 permFile <<
" " << checkPerm[i] << std::endl;
395 std::cout <<
"Checking permutation" << std::endl;
399 if (testReturn)
return testReturn;
403 std::cout <<
"Checking inverse permutation" << std::endl;
405 for (
size_t i=0; i< checkLength; i++){
406 testReturn = (checkInvPerm[checkPerm[i]] !=
z2TestLO(i));
407 if (testReturn)
return testReturn;
412 std::cout <<
"Checking num blocks" << std::endl;
414 testReturn = !((NumBlocks > 0) && (NumBlocks<
z2TestLO(checkLength)));
415 if (testReturn)
return testReturn;
420 std::cout <<
"Checking range" << std::endl;
422 testReturn = RangeTab[0];
423 if (testReturn)
return testReturn;
425 for (
z2TestLO i = 0; i < NumBlocks; i++){
426 testReturn = !(RangeTab[i] < RangeTab[i+1]);
427 if (testReturn)
return testReturn;
433 std::cout <<
"Checking Separator Tree" << std::endl;
437 testReturn = (TreeTab[0] != -1);
441 std::cout <<
"TreeTab[0] = " << TreeTab[0] <<
" != -1" << std::endl;
445 for (
size_t i=1; i< checkLength; i++){
446 testReturn = !(TreeTab[i] < NumBlocks);
448 std::cout <<
"TreeTab[" << i <<
"] = " << TreeTab[i] <<
" >= NumBlocks "
449 <<
" = " << NumBlocks << std::endl;
455 std::cout <<
"Going to compute the bandwidth" << std::endl;
463 std::cout <<
"Original Bandwidth: " << originalBandwidth << std::endl;
464 std::cout <<
"Permuted Bandwidth: " << permutedBandwidth << std::endl;
467 if(permutedBandwidth >= originalBandwidth) {
469 std::cout <<
"Test failed: bandwidth was not improved!" << std::endl;
471 std::cout <<
"Test in progress - returning OK for now..." << std::endl;
478 std::cout <<
"The bandwidth was improved. Good!" << std::endl;
482 catch (std::exception &e) {
483 std::cout <<
"Exception caught in ordering" << std::endl;
484 std::cout << e.what() << std::endl;
493 Tpetra::ScopeGuard tscope(&narg, &arg);
494 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
501 reduceAll<int,int>(*comm, Teuchos::EReductionType::REDUCE_MAX, 1,
502 &result, &resultReduced);
506 if (comm->getRank() == 0) {
507 std::cout <<
"Scotch Ordering test with " << comm->getSize()
508 <<
" processes has test return code " << resultReduced
509 <<
" and is exiting in the "
510 << ((resultReduced == 0 ) ?
"PASSED" :
"FAILED") <<
" state."
Provides access for Zoltan2 to Xpetra::CrsMatrix data.
Tpetra::CrsMatrix< z2TestScalar, z2TestLO, z2TestGO > SparseMatrix
int main(int narg, char **arg)
common code used by tests
size_t computeBandwidth(RCP< SparseMatrix > A, z2TestLO *iperm)
Tpetra::Vector< z2TestScalar, z2TestLO, z2TestGO > Vector
OrderingProblem sets up ordering problems for the user.
Defines the XpetraCrsMatrixAdapter class.
size_t getPermutationSize() const
Get (local) size of permutation.
lno_t getNumSeparatorBlocks() const
Get number of separator column blocks.
lno_t * getPermutationView(bool inverse=false) const
Get pointer to permutation. If inverse = true, return inverse permutation. By default, perm[i] is where new index i can be found in the old ordering. When inverse==true, perm[i] is where old index i can be found in the new ordering.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
map_t::local_ordinal_type lno_t
bool haveSeparators() const
Do we have the separators?
lno_t * getSeparatorRangeView() const
Get pointer to separator range.
Tpetra::Map::local_ordinal_type zlno_t
Defines the OrderingProblem class.
int validatePerm()
returns 0 if permutation is valid, negative if invalid.
Zoltan2::XpetraCrsMatrixAdapter< SparseMatrix > SparseMatrixAdapter
lno_t * getSeparatorTreeView() const
Get pointer to separator tree.
void tempDebugTest(RCP< SparseMatrix > origMatrix, Zoltan2::LocalOrderingSolution< z2TestLO > *soln)
bool mainExecute(int narg, char *arg[], RCP< const Comm< int > > &comm)
Tpetra::Map::global_ordinal_type zgno_t
typename InputTraits< User >::lno_t lno_t
LocalOrderingSolution< lno_t > * getLocalOrderingSolution()
Get the local ordering solution to the problem.
std::string testDataFilePath(".")