51 #include <Teuchos_ParameterList.hpp>
52 #include <Teuchos_RCP.hpp>
53 #include <Teuchos_CommandLineProcessor.hpp>
54 #include <Tpetra_CrsMatrix.hpp>
55 #include <Tpetra_Vector.hpp>
56 #include <MatrixMarket_Tpetra.hpp>
80 typedef Tpetra::CrsMatrix<z2TestScalar, z2TestLO, z2TestGO>
SparseMatrix;
81 typedef Tpetra::Vector<z2TestScalar, z2TestLO, z2TestGO>
Vector;
82 typedef Vector::node_type
Node;
86 #define epsilon 0.00000001
91 std::vector<int> count(n);
94 for (
size_t i=0; i<n; i++)
97 for (
size_t i=0; i<n; i++){
98 if (perm[i] < 0 || static_cast<size_t> (perm[i]) >= n)
105 for (
size_t i=0; i<n; i++){
118 std::string inputFile =
"";
119 std::string outputFile =
"";
120 bool verbose =
false;
124 Tpetra::ScopeGuard tscope(&narg, &arg);
125 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
126 int me = comm->getRank();
129 Teuchos::CommandLineProcessor cmdp (
false,
false);
130 cmdp.setOption(
"inputFile", &inputFile,
131 "Name of a Matrix Market file in the data directory; "
132 "if not specified, a matrix will be generated by MueLu.");
133 cmdp.setOption(
"outputFile", &outputFile,
134 "Name of the Matrix Market sparse matrix file to write, "
135 "echoing the input/generated matrix.");
136 cmdp.setOption(
"verbose",
"quiet", &verbose,
137 "Print messages and results.");
148 std::string matrixType(
"Laplace3D");
150 cmdp.setOption(
"x", &xdim,
151 "number of gridpoints in X dimension for "
152 "mesh used to generate matrix.");
153 cmdp.setOption(
"y", &ydim,
154 "number of gridpoints in Y dimension for "
155 "mesh used to generate matrix.");
156 cmdp.setOption(
"z", &zdim,
157 "number of gridpoints in Z dimension for "
158 "mesh used to generate matrix.");
159 cmdp.setOption(
"matrix", &matrixType,
160 "Matrix type: Laplace1D, Laplace2D, or Laplace3D");
163 cmdp.parse(narg, arg);
166 RCP<UserInputForTests> uinput;
172 uinput = rcp(
new UserInputForTests(xdim, ydim, zdim, matrixType, comm,
true,
true));
174 RCP<SparseMatrix> origMatrix = uinput->getUITpetraCrsMatrix();
176 if (outputFile !=
"") {
178 Tpetra::MatrixMarket::Writer<SparseMatrix>::writeSparseFile(outputFile,
179 origMatrix, verbose);
183 std::cout <<
"NumRows = " << origMatrix->getGlobalNumRows() << std::endl
184 <<
"NumNonzeros = " << origMatrix->getGlobalNumEntries() << std::endl
185 <<
"NumProcs = " << comm->getSize() << std::endl;
188 RCP<Vector> origVector, origProd;
189 origProd = Tpetra::createVector<z2TestScalar,z2TestLO,z2TestGO>(
190 origMatrix->getRangeMap());
191 origVector = Tpetra::createVector<z2TestScalar,z2TestLO,z2TestGO>(
192 origMatrix->getDomainMap());
193 origVector->randomize();
196 Teuchos::ParameterList params;
204 params.set(
"order_method",
"metis");
205 params.set(
"order_method_type",
"local");
220 for (
size_t ii = 0; ii < checkLength; ii++)
221 std::cout << checkPerm[ii] <<
" ";
222 std::cout << std::endl;
226 }
catch (std::exception &e){
227 #ifdef HAVE_ZOLTAN2_METIS
229 std::cout <<
"Exception from METIS Algorithm" << std::endl;
230 std::cout <<
"FAIL" << std::endl;
233 std::cout <<
"METIS is not enabled. METIS Algorithm threw an exception."
235 std::cout <<
"PASS" << std::endl;
242 std::cout <<
"Solution is not a permutation; FAIL" << std::endl;
244 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
Tpetra::Vector< z2TestScalar, z2TestLO, z2TestGO > Vector
OrderingProblem sets up ordering problems for the user.
int validatePerm(size_t n, z2TestLO *perm)
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.
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(".")