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;
90 typename SparseMatrix::local_inds_host_view_type indices;
91 typename SparseMatrix::values_host_view_type values;
99 for (ii=0; ii<n; ii++) {
100 A->getLocalRowView (ii, indices, values);
101 for (k=0; k< static_cast<z2TestLO>(indices.size()); k++){
105 j = iperm[indices[k]];
119 return (bw_left + bw_right + 1);
125 std::string inputFile =
"";
126 std::string outputFile =
"";
127 bool verbose =
false;
129 size_t origBandwidth = 0;
130 size_t permutedBandwidth = 0;
133 Tpetra::ScopeGuard tscope(&narg, &arg);
134 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
136 int me = comm->getRank();
139 Teuchos::CommandLineProcessor cmdp (
false,
false);
140 cmdp.setOption(
"inputFile", &inputFile,
141 "Name of a Matrix Market file in the data directory; "
142 "if not specified, a matrix will be generated by Galeri.");
143 cmdp.setOption(
"outputFile", &outputFile,
144 "Name of file to write the permutation");
145 cmdp.setOption(
"verbose",
"quiet", &verbose,
146 "Print messages and results.");
147 if (me == 0) std::cout <<
"Starting everything" << std::endl;
158 std::string matrixType(
"Laplace3D");
160 cmdp.setOption(
"x", &xdim,
161 "number of gridpoints in X dimension for "
162 "mesh used to generate matrix.");
163 cmdp.setOption(
"y", &ydim,
164 "number of gridpoints in Y dimension for "
165 "mesh used to generate matrix.");
166 cmdp.setOption(
"z", &zdim,
167 "number of gridpoints in Z dimension for "
168 "mesh used to generate matrix.");
169 cmdp.setOption(
"matrix", &matrixType,
170 "Matrix type: Laplace1D, Laplace2D, or Laplace3D");
175 std::string orderMethod(
"rcm");
176 cmdp.setOption(
"order_method", &orderMethod,
177 "order_method: natural, random, rcm, sorted_degree");
179 std::string orderMethodType(
"local");
180 cmdp.setOption(
"order_method_type", &orderMethodType,
181 "local or global or both");
184 cmdp.parse(narg, arg);
187 RCP<UserInputForTests> uinput;
189 if (inputFile !=
""){
194 uinput = rcp(
new UserInputForTests(xdim, ydim, zdim, matrixType, comm,
true,
true));
196 RCP<SparseMatrix> origMatrix = uinput->getUITpetraCrsMatrix();
199 std::cout <<
"NumRows = " << origMatrix->getGlobalNumRows() << std::endl
200 <<
"NumNonzeros = " << origMatrix->getGlobalNumEntries() << std::endl
201 <<
"NumProcs = " << comm->getSize() << std::endl;
215 Teuchos::ParameterList params;
216 params.set(
"order_method", orderMethod);
217 params.set(
"order_method_type", orderMethodType);
226 if (me == 0) std::cout <<
"Going to solve" << std::endl;
234 if (me == 0) std::cout <<
"Going to get results" << std::endl;
239 if (outputFile !=
"") {
240 std::ofstream permFile;
245 std::stringstream
fname;
246 fname << outputFile <<
"." << comm->getSize() <<
"." << me;
247 permFile.open(fname.str().c_str());
249 for (
size_t i=0; i<checkLength; i++){
250 permFile <<
" " << perm[i] << std::endl;
256 if (me == 0) std::cout <<
"Verifying results " << std::endl;
268 }
catch (std::exception &e){
269 std::cout <<
"Exception caught in ordering" << std::endl;
270 std::cout << e.what() << std::endl;
271 std::cout <<
"FAIL" << std::endl;
276 std::cout << me <<
": Not a valid permutation" << std::endl;
279 Teuchos::reduceAll<int,int>(*comm, Teuchos::EReductionType::REDUCE_MAX, 1,
280 &testReturn, &gTestReturn);
282 int increasedBandwidth = (permutedBandwidth > origBandwidth);
283 if (increasedBandwidth)
284 std::cout << me <<
": Bandwidth increased: original "
285 << origBandwidth <<
" < " << permutedBandwidth <<
" permuted "
288 std::cout << me <<
": Bandwidth not increased: original "
289 << origBandwidth <<
" >= " << permutedBandwidth <<
" permuted "
292 int gIncreasedBandwidth;
293 Teuchos::reduceAll<int,int>(*comm, Teuchos::EReductionType::REDUCE_MAX, 1,
294 &increasedBandwidth, &gIncreasedBandwidth);
298 std::cout <<
"Solution is not a permutation; FAIL" << std::endl;
299 else if (gIncreasedBandwidth && (orderMethod ==
"rcm"))
300 std::cout <<
"Bandwidth was increased; FAIL" << std::endl;
302 std::cout <<
"PASS" << std::endl;
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 * 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.
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
Tpetra::Map::global_ordinal_type zgno_t
LocalOrderingSolution< lno_t > * getLocalOrderingSolution()
Get the local ordering solution to the problem.
std::string testDataFilePath(".")