Zoltan2
XpetraCrsMatrixInput.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
45 //
46 // Basic testing of Zoltan2::XpetraCrsMatrixAdapter
47 
53 #include <string>
54 
56 #include <Zoltan2_InputTraits.hpp>
57 #include <Zoltan2_TestHelpers.hpp>
58 
59 #include <Teuchos_GlobalMPISession.hpp>
60 #include <Teuchos_DefaultComm.hpp>
61 #include <Teuchos_RCP.hpp>
62 #include <Teuchos_Comm.hpp>
63 #include <Teuchos_CommHelpers.hpp>
64 
65 using namespace std;
66 using Teuchos::RCP;
67 using Teuchos::rcp;
68 using Teuchos::rcp_const_cast;
69 using Teuchos::Comm;
70 using Teuchos::DefaultComm;
71 
72 typedef Tpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t> tmatrix_t;
73 typedef Xpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t> xmatrix_t;
74 
75 void printMatrix(RCP<const Comm<int> > &comm, zlno_t nrows,
76  const zgno_t *rowIds, const zlno_t *offsets, const zgno_t *colIds)
77 {
78  int rank = comm->getRank();
79  int nprocs = comm->getSize();
80  comm->barrier();
81  for (int p=0; p < nprocs; p++){
82  if (p == rank){
83  std::cout << rank << ":" << std::endl;
84  for (zlno_t i=0; i < nrows; i++){
85  std::cout << " row " << rowIds[i] << ": ";
86  for (zlno_t j=offsets[i]; j < offsets[i+1]; j++){
87  std::cout << colIds[j] << " ";
88  }
89  std::cout << std::endl;
90  }
91  std::cout.flush();
92  }
93  comm->barrier();
94  }
95  comm->barrier();
96 }
97 
98 template <typename User>
101 {
102  RCP<const Comm<int> > comm = M.getComm();
103  int fail = 0, gfail=0;
104 
105  if (!fail && ia.getLocalNumRows() != M.getNodeNumRows())
106  fail = 4;
107 
108  if (M.getNodeNumRows()){
109  if (!fail && ia.getLocalNumColumns() != M.getNodeNumCols())
110  fail = 6;
111  }
112 
113  gfail = globalFail(comm, fail);
114 
115  const zgno_t *rowIds=NULL, *colIds=NULL;
116  const zlno_t *offsets=NULL;
117  size_t nrows=0;
118 
119  if (!gfail){
120 
121  nrows = ia.getLocalNumRows();
122  ia.getRowIDsView(rowIds);
123  ia.getCRSView(offsets, colIds);
124 
125  if (nrows != M.getNodeNumRows())
126  fail = 8;
127 
128  gfail = globalFail(comm, fail);
129 
130  if (gfail == 0){
131  printMatrix(comm, nrows, rowIds, offsets, colIds);
132  }
133  else{
134  if (!fail) fail = 10;
135  }
136  }
137  return fail;
138 }
139 
140 int main(int argc, char *argv[])
141 {
142  Teuchos::GlobalMPISession session(&argc, &argv);
143  RCP<const Comm<int> > comm = DefaultComm<int>::getComm();
144  int rank = comm->getRank();
145  int fail = 0, gfail=0;
146  bool aok = true;
147 
148  // Create object that can give us test Tpetra, Xpetra
149  // and Epetra matrices for testing.
150 
151  RCP<UserInputForTests> uinput;
152 
153  try{
154  uinput =
155  rcp(new UserInputForTests(
156  testDataFilePath,std::string("simple"), comm, true));
157  }
158  catch(std::exception &e){
159  aok = false;
160  std::cout << e.what() << std::endl;
161  }
162  TEST_FAIL_AND_EXIT(*comm, aok, "input ", 1);
163 
164  RCP<tmatrix_t> tM; // original matrix (for checking)
165  RCP<tmatrix_t> newM; // migrated matrix
166 
167  tM = uinput->getUITpetraCrsMatrix();
168  size_t nrows = tM->getNodeNumRows();
169 
170  // To test migration in the input adapter we need a Solution object.
171 
172  RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment);
173 
174  int nWeights = 1;
175 
178  typedef adapter_t::part_t part_t;
179 
180  part_t *p = new part_t [nrows];
181  memset(p, 0, sizeof(part_t) * nrows);
182  ArrayRCP<part_t> solnParts(p, 0, nrows, true);
183 
184  soln_t solution(env, comm, nWeights);
185  solution.setParts(solnParts);
186 
188  // User object is Tpetra::CrsMatrix
189  if (!gfail){
190  RCP<const tmatrix_t> ctM = rcp_const_cast<const tmatrix_t>(tM);
191  RCP<Zoltan2::XpetraCrsMatrixAdapter<tmatrix_t> > tMInput;
192 
193  try {
194  tMInput =
196  }
197  catch (std::exception &e){
198  aok = false;
199  std::cout << e.what() << std::endl;
200  }
201  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter ", 1);
202 
203  if (rank==0)
204  std::cout << "Input adapter for Tpetra::CrsMatrix" << std::endl;
205 
206  fail = verifyInputAdapter<tmatrix_t>(*tMInput, *tM);
207 
208  gfail = globalFail(comm, fail);
209 
210  if (!gfail){
211  tmatrix_t *mMigrate = NULL;
212  try{
213  tMInput->applyPartitioningSolution(*tM, mMigrate, solution);
214  newM = rcp(mMigrate);
215  }
216  catch (std::exception &e){
217  fail = 11;
218  cout << "Error caught: " << e.what() << endl;
219  }
220 
221  gfail = globalFail(comm, fail);
222 
223  if (!gfail){
224  RCP<const tmatrix_t> cnewM = rcp_const_cast<const tmatrix_t>(newM);
225  RCP<Zoltan2::XpetraCrsMatrixAdapter<tmatrix_t> > newInput;
226  try{
227  newInput = rcp(new Zoltan2::XpetraCrsMatrixAdapter<tmatrix_t>(cnewM));
228  }
229  catch (std::exception &e){
230  aok = false;
231  std::cout << e.what() << std::endl;
232  }
233  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter 2 ", 1);
234 
235  if (rank==0){
236  std::cout <<
237  "Input adapter for Tpetra::CrsMatrix migrated to proc 0" <<
238  std::endl;
239  }
240  fail = verifyInputAdapter<tmatrix_t>(*newInput, *newM);
241  if (fail) fail += 100;
242  gfail = globalFail(comm, fail);
243  }
244  }
245  if (gfail){
246  printFailureCode(comm, fail);
247  }
248  }
249 
251  // User object is Xpetra::CrsMatrix
252  if (!gfail){
253  RCP<xmatrix_t> xM = uinput->getUIXpetraCrsMatrix();
254  RCP<const xmatrix_t> cxM = rcp_const_cast<const xmatrix_t>(xM);
255  RCP<Zoltan2::XpetraCrsMatrixAdapter<xmatrix_t> > xMInput;
256 
257  try {
258  xMInput =
260  }
261  catch (std::exception &e){
262  aok = false;
263  std::cout << e.what() << std::endl;
264  }
265  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter 3 ", 1);
266 
267  if (rank==0){
268  std::cout << "Input adapter for Xpetra::CrsMatrix" << std::endl;
269  }
270  fail = verifyInputAdapter<xmatrix_t>(*xMInput, *tM);
271 
272  gfail = globalFail(comm, fail);
273 
274  if (!gfail){
275  xmatrix_t *mMigrate =NULL;
276  try{
277  xMInput->applyPartitioningSolution(*xM, mMigrate, solution);
278  }
279  catch (std::exception &e){
280  cout << "Error caught: " << e.what() << endl;
281  fail = 11;
282  }
283 
284  gfail = globalFail(comm, fail);
285 
286  if (!gfail){
287  RCP<const xmatrix_t> cnewM(mMigrate);
288  RCP<Zoltan2::XpetraCrsMatrixAdapter<xmatrix_t> > newInput;
289  try{
290  newInput =
292  }
293  catch (std::exception &e){
294  aok = false;
295  std::cout << e.what() << std::endl;
296  }
297  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter 4 ", 1);
298 
299  if (rank==0){
300  std::cout <<
301  "Input adapter for Xpetra::CrsMatrix migrated to proc 0" <<
302  std::endl;
303  }
304  fail = verifyInputAdapter<xmatrix_t>(*newInput, *newM);
305  if (fail) fail += 100;
306  gfail = globalFail(comm, fail);
307  }
308  }
309  if (gfail){
310  printFailureCode(comm, fail);
311  }
312  }
313 
314 #ifdef HAVE_EPETRA_DATA_TYPES
315  // User object is Epetra_CrsMatrix
317  typedef Epetra_CrsMatrix ematrix_t;
318  if (!gfail){
319  RCP<ematrix_t> eM = uinput->getUIEpetraCrsMatrix();
320  RCP<const ematrix_t> ceM = rcp_const_cast<const ematrix_t>(eM);
321  RCP<Zoltan2::XpetraCrsMatrixAdapter<ematrix_t> > eMInput;
322 
323  try {
324  eMInput =
326  }
327  catch (std::exception &e){
328  aok = false;
329  std::cout << e.what() << std::endl;
330  }
331  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter 5 ", 1);
332 
333  if (rank==0){
334  std::cout << "Input adapter for Epetra_CrsMatrix" << std::endl;
335  }
336  fail = verifyInputAdapter<ematrix_t>(*eMInput, *tM);
337 
338  gfail = globalFail(comm, fail);
339 
340  if (!gfail){
341  ematrix_t *mMigrate =NULL;
342  try{
343  eMInput->applyPartitioningSolution(*eM, mMigrate, solution);
344  }
345  catch (std::exception &e){
346  cout << "Error caught: " << e.what() << endl;
347  fail = 11;
348  }
349 
350  gfail = globalFail(comm, fail);
351 
352  if (!gfail){
353  RCP<const ematrix_t> cnewM(mMigrate, true);
354  RCP<Zoltan2::XpetraCrsMatrixAdapter<ematrix_t> > newInput;
355  try{
356  newInput =
358  }
359  catch (std::exception &e){
360  aok = false;
361  std::cout << e.what() << std::endl;
362  }
363  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter 6 ", 1);
364 
365  if (rank==0){
366  std::cout <<
367  "Input adapter for Epetra_CrsMatrix migrated to proc 0" <<
368  std::endl;
369  }
370  fail = verifyInputAdapter<ematrix_t>(*newInput, *newM);
371  if (fail) fail += 100;
372  gfail = globalFail(comm, fail);
373  }
374  }
375  if (gfail){
376  printFailureCode(comm, fail);
377  }
378  }
379 #endif
380 
382  // DONE
383 
384  if (rank==0)
385  std::cout << "PASS" << std::endl;
386 }
int globalFail(const RCP< const Comm< int > > &comm, int fail)
Provides access for Zoltan2 to Xpetra::CrsMatrix data.
Tpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > tmatrix_t
Xpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > xmatrix_t
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
void getCRSView(const lno_t *&offsets, const gno_t *&colIds) const
int main(int argc, char *argv[])
size_t getLocalNumColumns() const
Returns the number of columns on this process.
int zlno_t
common code used by tests
Defines the XpetraCrsMatrixAdapter class.
void printMatrix(RCP< const Comm< int > > &comm, zlno_t nrows, const zgno_t *rowIds, const zlno_t *offsets, const zgno_t *colIds)
A PartitioningSolution is a solution to a partitioning problem.
size_t getLocalNumRows() const
Returns the number of rows on this process.
Traits for application input objects.
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
static const std::string fail
int zgno_t
int verifyInputAdapter(Zoltan2::XpetraCrsMatrixAdapter< User > &ia, tmatrix_t &M)
void printFailureCode(const RCP< const Comm< int > > &comm, int fail)
void getRowIDsView(const gno_t *&rowIds) const
std::string testDataFilePath(".")