Zoltan2
Loading...
Searching...
No Matches
ordering1.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
48#include <iostream>
49#include <fstream>
50#include <limits>
51#include <vector>
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>
58
59using Teuchos::RCP;
60
62// Program to demonstrate use of Zoltan2 to order a TPetra matrix
63// (read from a MatrixMarket file or generated by Galeri::Xpetra).
64// Usage:
65// a.out [--inputFile=filename] [--outputFile=outfile] [--verbose]
66// [--x=#] [--y=#] [--z=#] [--matrix={Laplace1D,Laplace2D,Laplace3D}
67// Karen Devine, 2011
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;
83
84size_t computeBandwidth(RCP<SparseMatrix> A, z2TestLO *iperm)
85// Returns the bandwidth of the (local) permuted matrix
86// Uses the inverse permutation calculated from the OrderingSolution
87// if passed in, otherwise is calculating the original value.
88{
89 z2TestLO ii, i, j, k;
90 typename SparseMatrix::local_inds_host_view_type indices;
91 typename SparseMatrix::values_host_view_type values;
92
93 z2TestLO bw_left = 0;
94 z2TestLO bw_right = 0;
95
96 z2TestLO n = A->getLocalNumRows();
97
98 // Loop over rows of matrix
99 for (ii=0; ii<n; ii++) {
100 A->getLocalRowView (ii, indices, values);
101 for (k=0; k< static_cast<z2TestLO>(indices.size()); k++){
102 if (indices[k] < n){ // locally owned
103 if (iperm){
104 i = iperm[ii];
105 j = iperm[indices[k]];
106 } else {
107 i = ii;
108 j = indices[k];
109 }
110 if (j-i > bw_right)
111 bw_right = j-i;
112 if (i-j > bw_left)
113 bw_left = i-j;
114 }
115 }
116 }
117
118 // Total bandwidth is the sum of left and right + 1
119 return (bw_left + bw_right + 1);
120}
121
123int main(int narg, char** arg)
124{
125 std::string inputFile = ""; // Matrix Market file to read
126 std::string outputFile = ""; // Output file to write
127 bool verbose = false; // Verbosity of output
128 int testReturn = 0;
129 size_t origBandwidth = 0;
130 size_t permutedBandwidth = 0;
131
133 Tpetra::ScopeGuard tscope(&narg, &arg);
134 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
135
136 int me = comm->getRank();
137
138 // Read run-time options.
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;
148
150 // Even with cmdp option "true", I get errors for having these
151 // arguments on the command line. (On redsky build)
152 // KDDKDD Should just be warnings, right? Code should still work with these
153 // KDDKDD params in the create-a-matrix file. Better to have them where
154 // KDDKDD they are used.
155 int xdim=10;
156 int ydim=10;
157 int zdim=10;
158 std::string matrixType("Laplace3D");
159
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");
171
173 // Ordering options to test.
175 std::string orderMethod("rcm"); // TODO: Allow "RCM" as well
176 cmdp.setOption("order_method", &orderMethod,
177 "order_method: natural, random, rcm, sorted_degree");
178
179 std::string orderMethodType("local"); // TODO: Allow "LOCAL" as well
180 cmdp.setOption("order_method_type", &orderMethodType,
181 "local or global or both");
182
184 cmdp.parse(narg, arg);
185
186
187 RCP<UserInputForTests> uinput;
188
189 if (inputFile != ""){ // Input file specified; read a matrix
190 uinput = rcp(new UserInputForTests(testDataFilePath, inputFile, comm, true));
191 }
192 else // Let Galeri generate a matrix
193
194 uinput = rcp(new UserInputForTests(xdim, ydim, zdim, matrixType, comm, true, true));
195
196 RCP<SparseMatrix> origMatrix = uinput->getUITpetraCrsMatrix();
197
198 if (me == 0)
199 std::cout << "NumRows = " << origMatrix->getGlobalNumRows() << std::endl
200 << "NumNonzeros = " << origMatrix->getGlobalNumEntries() << std::endl
201 << "NumProcs = " << comm->getSize() << std::endl;
202
204 // Currently Not Used
205 /*
206 RCP<Vector> origVector, origProd;
207 origProd = Tpetra::createVector<z2TestScalar,z2TestLO,z2TestGO>(
208 origMatrix->getRangeMap());
209 origVector = Tpetra::createVector<z2TestScalar,z2TestLO,z2TestGO>(
210 origMatrix->getDomainMap());
211 origVector->randomize();
212 */
213
215 Teuchos::ParameterList params;
216 params.set("order_method", orderMethod);
217 params.set("order_method_type", orderMethodType);
218
220 SparseMatrixAdapter adapter(origMatrix);
221
223 try
224 {
225 Zoltan2::OrderingProblem<SparseMatrixAdapter> problem(&adapter, &params);
226 if (me == 0) std::cout << "Going to solve" << std::endl;
227 problem.solve();
228
230
232 problem.getLocalOrderingSolution();
233
234 if (me == 0) std::cout << "Going to get results" << std::endl;
235 // Check that the solution is really a permutation
236
237 z2TestLO * perm = soln->getPermutationView();
238
239 if (outputFile != "") {
240 std::ofstream permFile;
241
242 // Write permutation (0-based) to file
243 // each process writes local perm to a separate file
244 //std::string fname = outputFile + "." + me;
245 std::stringstream fname;
246 fname << outputFile << "." << comm->getSize() << "." << me;
247 permFile.open(fname.str().c_str());
248 size_t checkLength = soln->getPermutationSize();
249 for (size_t i=0; i<checkLength; i++){
250 permFile << " " << perm[i] << std::endl;
251 }
252 permFile.close();
253
254 }
255
256 if (me == 0) std::cout << "Verifying results " << std::endl;
257
258 // Verify that checkPerm is a permutation
259 testReturn = soln->validatePerm();
260
261 // Compute original bandwidth
262 origBandwidth = computeBandwidth(origMatrix, nullptr);
263
264 // Compute permuted bandwidth
265 z2TestLO * iperm = soln->getPermutationView(true);
266 permutedBandwidth = computeBandwidth(origMatrix, iperm);
267
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;
272 return 0;
273 }
274
275 if (testReturn)
276 std::cout << me << ": Not a valid permutation" << std::endl;
277
278 int gTestReturn;
279 Teuchos::reduceAll<int,int>(*comm, Teuchos::EReductionType::REDUCE_MAX, 1,
280 &testReturn, &gTestReturn);
281
282 int increasedBandwidth = (permutedBandwidth > origBandwidth);
283 if (increasedBandwidth)
284 std::cout << me << ": Bandwidth increased: original "
285 << origBandwidth << " < " << permutedBandwidth << " permuted "
286 << std::endl;
287 else
288 std::cout << me << ": Bandwidth not increased: original "
289 << origBandwidth << " >= " << permutedBandwidth << " permuted "
290 << std::endl;
291
292 int gIncreasedBandwidth;
293 Teuchos::reduceAll<int,int>(*comm, Teuchos::EReductionType::REDUCE_MAX, 1,
294 &increasedBandwidth, &gIncreasedBandwidth);
295
296 if (me == 0) {
297 if (gTestReturn)
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;
301 else
302 std::cout << "PASS" << std::endl;
303
304 }
305
306}
307
Defines the OrderingProblem 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()
OrderingProblem sets up ordering problems for the user.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
LocalOrderingSolution< lno_t > * getLocalOrderingSolution()
Get the local ordering solution to the problem.
int validatePerm()
returns 0 if permutation is valid, negative if invalid.
index_t * getPermutationView(bool inverse=false) const
Get pointer to permutation. If inverse = true, return inverse permutation. By default,...
size_t getPermutationSize() const
Get (local) size of permutation.
Provides access for Zoltan2 to Xpetra::CrsMatrix data.
zlno_t z2TestLO
Definition: coloring1.cpp:75
size_t computeBandwidth(RCP< SparseMatrix > A, z2TestLO *iperm)
Definition: ordering1.cpp:84
Vector::node_type Node
Definition: ordering1.cpp:81
zlno_t z2TestLO
Definition: ordering1.cpp:75
Tpetra::CrsMatrix< z2TestScalar, z2TestLO, z2TestGO > SparseMatrix
Definition: ordering1.cpp:79
Zoltan2::XpetraCrsMatrixAdapter< SparseMatrix > SparseMatrixAdapter
Definition: ordering1.cpp:82
Tpetra::Vector< z2TestScalar, z2TestLO, z2TestGO > Vector
Definition: ordering1.cpp:80
zgno_t z2TestGO
Definition: ordering1.cpp:76
zscalar_t z2TestScalar
Definition: ordering1.cpp:77