Zoltan2
XpetraCrsGraphInput.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::XpetraCrsGraphAdapter
52 #include <string>
53 
56 #include <Zoltan2_TestHelpers.hpp>
57 
58 #include <Teuchos_GlobalMPISession.hpp>
59 #include <Teuchos_DefaultComm.hpp>
60 #include <Teuchos_RCP.hpp>
61 #include <Teuchos_Comm.hpp>
62 #include <Teuchos_CommHelpers.hpp>
63 
64 using namespace std;
65 using Teuchos::RCP;
66 using Teuchos::rcp;
67 using Teuchos::rcp_const_cast;
68 using Teuchos::Comm;
69 using Teuchos::DefaultComm;
70 using Teuchos::Array;
71 using Teuchos::ArrayView;
72 
73 typedef Tpetra::CrsGraph<zlno_t, zgno_t, znode_t> tgraph_t;
74 typedef Xpetra::CrsGraph<zlno_t, zgno_t, znode_t> xgraph_t;
75 
76 void printGraph(RCP<const Comm<int> > &comm, zlno_t nvtx,
77  const zgno_t *vtxIds, const zlno_t *offsets, const zgno_t *edgeIds)
78 {
79  int rank = comm->getRank();
80  int nprocs = comm->getSize();
81  comm->barrier();
82  for (int p=0; p < nprocs; p++){
83  if (p == rank){
84  std::cout << rank << ":" << std::endl;
85  for (zlno_t i=0; i < nvtx; i++){
86  std::cout << " vertex " << vtxIds[i] << ": ";
87  for (zlno_t j=offsets[i]; j < offsets[i+1]; j++){
88  std::cout << edgeIds[j] << " ";
89  }
90  std::cout << std::endl;
91  }
92  std::cout.flush();
93  }
94  comm->barrier();
95  }
96  comm->barrier();
97 }
98 
99 template <typename User>
102  tgraph_t &graph
103 )
104 {
105  RCP<const Comm<int> > comm = graph.getComm();
106  int fail = 0, gfail=0;
107 
108  if (!fail &&
109  ia.getLocalNumVertices() != graph.getNodeNumRows())
110  fail = 4;
111 
112  if (!fail &&
113  ia.getLocalNumEdges() != graph.getNodeNumEntries())
114  fail = 6;
115 
116  gfail = globalFail(comm, fail);
117 
118  const zgno_t *vtxIds=NULL, *edgeIds=NULL;
119  const zlno_t *offsets=NULL;
120  size_t nvtx=0;
121 
122  if (!gfail){
123 
124  nvtx = ia.getLocalNumVertices();
125  ia.getVertexIDsView(vtxIds);
126  ia.getEdgesView(offsets, edgeIds);
127 
128  if (nvtx != graph.getNodeNumRows())
129  fail = 8;
130 
131  gfail = globalFail(comm, fail);
132 
133  if (gfail == 0){
134  printGraph(comm, nvtx, vtxIds, offsets, edgeIds);
135  }
136  else{
137  if (!fail) fail = 10;
138  }
139  }
140  return fail;
141 }
142 
143 int main(int argc, char *argv[])
144 {
145  Teuchos::GlobalMPISession session(&argc, &argv);
146  RCP<const Comm<int> > comm = DefaultComm<int>::getComm();
147  int rank = comm->getRank();
148  int fail = 0, gfail=0;
149  bool aok = true;
150 
151  // Create an object that can give us test Tpetra, Xpetra
152  // and Epetra graphs for testing.
153 
154  RCP<UserInputForTests> uinput;
155 
156  try{
157  uinput =
158  rcp(new UserInputForTests(testDataFilePath,std::string("simple"), comm, true));
159  }
160  catch(std::exception &e){
161  aok = false;
162  std::cout << e.what() << std::endl;
163  }
164  TEST_FAIL_AND_EXIT(*comm, aok, "input ", 1);
165 
166  RCP<tgraph_t> tG; // original graph (for checking)
167  RCP<tgraph_t> newG; // migrated graph
168 
169  tG = uinput->getUITpetraCrsGraph();
170  size_t nvtx = tG->getNodeNumRows();
171 
172  // To test migration in the input adapter we need a Solution object.
173  // Our solution just assigns all objects to part zero.
174 
175  RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment);
176 
177  int nWeights = 1;
178 
181  typedef adapter_t::part_t part_t;
182 
183  part_t *p = new part_t [nvtx];
184  memset(p, 0, sizeof(part_t) * nvtx);
185  ArrayRCP<part_t> solnParts(p, 0, nvtx, true);
186 
187  soln_t solution(env, comm, nWeights);
188  solution.setParts(solnParts);
189 
191  // User object is Tpetra::CrsGraph
192  if (!gfail){
193  if (rank==0)
194  std::cout << "Input adapter for Tpetra::CrsGraph" << std::endl;
195 
196  RCP<const tgraph_t> ctG = rcp_const_cast<const tgraph_t>(tG);
197  RCP<Zoltan2::XpetraCrsGraphAdapter<tgraph_t> > tGInput;
198 
199  try {
200  tGInput =
202  }
203  catch (std::exception &e){
204  aok = false;
205  std::cout << e.what() << std::endl;
206  }
207  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsGraphAdapter ", 1);
208 
209  fail = verifyInputAdapter<tgraph_t>(*tGInput, *tG);
210 
211  gfail = globalFail(comm, fail);
212 
213  if (!gfail){
214  tgraph_t *mMigrate = NULL;
215  try{
216  tGInput->applyPartitioningSolution( *tG, mMigrate, solution);
217  newG = rcp(mMigrate);
218  }
219  catch (std::exception &e){
220  fail = 11;
221  }
222 
223  gfail = globalFail(comm, fail);
224 
225  if (!gfail){
226  RCP<const tgraph_t> cnewG = rcp_const_cast<const tgraph_t>(newG);
227  RCP<Zoltan2::XpetraCrsGraphAdapter<tgraph_t> > newInput;
228  try{
229  newInput = rcp(new Zoltan2::XpetraCrsGraphAdapter<tgraph_t>(cnewG));
230  }
231  catch (std::exception &e){
232  aok = false;
233  std::cout << e.what() << std::endl;
234  }
235  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsGraphAdapter 2 ", 1);
236 
237  if (rank==0){
238  std::cout <<
239  "Input adapter for Tpetra::CrsGraph migrated to proc 0" <<
240  std::endl;
241  }
242  fail = verifyInputAdapter<tgraph_t>(*newInput, *newG);
243  if (fail) fail += 100;
244  gfail = globalFail(comm, fail);
245  }
246  }
247  if (gfail){
248  printFailureCode(comm, fail);
249  }
250  }
251 
253  // User object is Xpetra::CrsGraph
254  if (!gfail){
255  if (rank==0)
256  std::cout << "Input adapter for Xpetra::CrsGraph" << std::endl;
257 
258  RCP<xgraph_t> xG = uinput->getUIXpetraCrsGraph();
259  RCP<const xgraph_t> cxG = rcp_const_cast<const xgraph_t>(xG);
260  RCP<Zoltan2::XpetraCrsGraphAdapter<xgraph_t> > xGInput;
261 
262  try {
263  xGInput =
265  }
266  catch (std::exception &e){
267  aok = false;
268  std::cout << e.what() << std::endl;
269  }
270  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsGraphAdapter 3 ", 1);
271 
272  fail = verifyInputAdapter<xgraph_t>(*xGInput, *tG);
273 
274  gfail = globalFail(comm, fail);
275 
276  if (!gfail){
277  xgraph_t *mMigrate =NULL;
278  try{
279  xGInput->applyPartitioningSolution( *xG, mMigrate, solution);
280  }
281  catch (std::exception &e){
282  fail = 11;
283  }
284 
285  gfail = globalFail(comm, fail);
286 
287  if (!gfail){
288  RCP<const xgraph_t> cnewG(mMigrate);
289  RCP<Zoltan2::XpetraCrsGraphAdapter<xgraph_t> > newInput;
290  try{
291  newInput =
293  }
294  catch (std::exception &e){
295  aok = false;
296  std::cout << e.what() << std::endl;
297  }
298  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsGraphAdapter 4 ", 1);
299 
300  if (rank==0){
301  std::cout <<
302  "Input adapter for Xpetra::CrsGraph migrated to proc 0" <<
303  std::endl;
304  }
305  fail = verifyInputAdapter<xgraph_t>(*newInput, *newG);
306  if (fail) fail += 100;
307  gfail = globalFail(comm, fail);
308  }
309  }
310  if (gfail){
311  printFailureCode(comm, fail);
312  }
313  }
314 
315 #ifdef HAVE_EPETRA_DATA_TYPES
316  // User object is Epetra_CrsGraph
318  typedef Epetra_CrsGraph egraph_t;
319  if (!gfail){
320  if (rank==0)
321  std::cout << "Input adapter for Epetra_CrsGraph" << std::endl;
322 
323  RCP<egraph_t> eG = uinput->getUIEpetraCrsGraph();
324  RCP<const egraph_t> ceG = rcp_const_cast<const egraph_t>(eG);
325  RCP<Zoltan2::XpetraCrsGraphAdapter<egraph_t> > eGInput;
326 
327  bool goodAdapter = true;
328  try {
329  eGInput =
331  }
332  catch (std::exception &e){
333  if (std::is_same<znode_t, Xpetra::EpetraNode>::value) {
334  aok = false;
335  goodAdapter = false;
336  std::cout << e.what() << std::endl;
337  }
338  else {
339  // We expect an error from Xpetra when znode_t != Xpetra::EpetraNode
340  // Ignore it, but skip tests using this matrix adapter.
341  std::cout << "Node type is not supported by Xpetra's Epetra interface;"
342  << " Skipping this test." << std::endl;
343  std::cout << "FYI: Here's the exception message: " << std::endl
344  << e.what() << std::endl;
345  goodAdapter = false;
346  }
347  }
348  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsGraphAdapter 5 ", 1);
349 
350  if (goodAdapter) {
351  fail = verifyInputAdapter<egraph_t>(*eGInput, *tG);
352 
353  gfail = globalFail(comm, fail);
354 
355  if (!gfail){
356  egraph_t *mMigrate =NULL;
357  try{
358  eGInput->applyPartitioningSolution( *eG, mMigrate, solution);
359  }
360  catch (std::exception &e){
361  fail = 11;
362  }
363 
364  gfail = globalFail(comm, fail);
365 
366  if (!gfail){
367  RCP<const egraph_t> cnewG(mMigrate, true);
368  RCP<Zoltan2::XpetraCrsGraphAdapter<egraph_t> > newInput;
369  try{
370  newInput =
372  }
373  catch (std::exception &e){
374  aok = false;
375  std::cout << e.what() << std::endl;
376  }
377  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsGraphAdapter 6 ", 1);
378 
379  if (rank==0){
380  std::cout <<
381  "Input adapter for Epetra_CrsGraph migrated to proc 0" <<
382  std::endl;
383  }
384  fail = verifyInputAdapter<egraph_t>(*newInput, *newG);
385  if (fail) fail += 100;
386  gfail = globalFail(comm, fail);
387  }
388  }
389  if (gfail){
390  printFailureCode(comm, fail);
391  }
392  }
393  }
394 #endif
395 
397  // DONE
398 
399  if (rank==0)
400  std::cout << "PASS" << std::endl;
401 }
void getVertexIDsView(const gno_t *&ids) const
Sets pointers to this process&#39; graph entries.
int globalFail(const RCP< const Comm< int > > &comm, int fail)
void printGraph(RCP< const Comm< int > > &comm, zlno_t nvtx, const zgno_t *vtxIds, const zlno_t *offsets, const zgno_t *edgeIds)
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
Provides access for Zoltan2 to Xpetra::CrsGraph data.
size_t getLocalNumEdges() const
Returns the number of edges on this process.
Xpetra::CrsGraph< zlno_t, zgno_t, znode_t > xgraph_t
int zlno_t
Defines the PartitioningSolution class.
common code used by tests
size_t getLocalNumVertices() const
Returns the number of vertices on this process.
Defines XpetraCrsGraphAdapter class.
A PartitioningSolution is a solution to a partitioning problem.
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::XpetraCrsGraphAdapter< User > &ia, tgraph_t &graph)
void printFailureCode(const RCP< const Comm< int > > &comm, int fail)
int main(int argc, char *argv[])
void getEdgesView(const lno_t *&offsets, const gno_t *&adjIds) const
Gets adjacency lists for all vertices in a compressed sparse row (CSR) format.
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > tgraph_t
std::string testDataFilePath(".")