Zoltan2
TpetraRowGraphInput.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::TpetraRowGraphAdapter
52 #include <string>
53 
55 #include <Zoltan2_InputTraits.hpp>
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::rcp_dynamic_cast;
69 using Teuchos::Comm;
70 using Teuchos::DefaultComm;
71 
72 typedef Tpetra::CrsGraph<zlno_t, zgno_t, znode_t> ztcrsgraph_t;
73 typedef Tpetra::RowGraph<zlno_t, zgno_t, znode_t> ztrowgraph_t;
74 
75 void printGraph(RCP<const Comm<int> > &comm, zlno_t nvtx,
76  const zgno_t *vtxIds, const zlno_t *offsets, const zgno_t *edgeIds)
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 < nvtx; i++){
85  std::cout << " vertex " << vtxIds[i] << ": ";
86  for (zlno_t j=offsets[i]; j < offsets[i+1]; j++){
87  std::cout << edgeIds[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 = graph.getComm();
103  int fail = 0, gfail=0;
104 
105  if (!fail &&
106  ia.getLocalNumVertices() != graph.getNodeNumRows())
107  fail = 4;
108 
109  if (!fail &&
110  ia.getLocalNumEdges() != graph.getNodeNumEntries())
111  fail = 6;
112 
113  gfail = globalFail(comm, fail);
114 
115  const zgno_t *vtxIds=NULL, *edgeIds=NULL;
116  const zlno_t *offsets=NULL;
117  size_t nvtx=0;
118 
119  if (!gfail){
120 
121  nvtx = ia.getLocalNumVertices();
122  ia.getVertexIDsView(vtxIds);
123  ia.getEdgesView(offsets, edgeIds);
124 
125  if (nvtx != graph.getNodeNumRows())
126  fail = 8;
127 
128  gfail = globalFail(comm, fail);
129 
130  if (gfail == 0){
131  printGraph(comm, nvtx, vtxIds, offsets, edgeIds);
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 an object that can give us test Tpetra graphs for testing
149 
150  RCP<UserInputForTests> uinput;
151 
152  try{
153  uinput =
154  rcp(new UserInputForTests(
155  testDataFilePath,std::string("simple"), comm, true));
156  }
157  catch(std::exception &e){
158  aok = false;
159  std::cout << e.what() << std::endl;
160  }
161  TEST_FAIL_AND_EXIT(*comm, aok, "input ", 1);
162 
163  // Input crs graph and row graph cast from it.
164  RCP<ztcrsgraph_t> tG = uinput->getUITpetraCrsGraph();
165  RCP<ztrowgraph_t> trG = rcp_dynamic_cast<ztrowgraph_t>(tG);
166 
167  RCP<ztrowgraph_t> newG; // migrated graph
168 
169  size_t nvtx = tG->getNodeNumRows();
170 
171  // To test migration in the input adapter we need a Solution object.
172 
173  RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment);
174 
175  int nWeights = 1;
176 
179  typedef adapter_t::part_t part_t;
180 
181  part_t *p = new part_t [nvtx];
182  memset(p, 0, sizeof(part_t) * nvtx);
183  ArrayRCP<part_t> solnParts(p, 0, nvtx, true);
184 
185  soln_t solution(env, comm, nWeights);
186  solution.setParts(solnParts);
187 
189  // User object is Tpetra::RowGraph
190  if (!gfail){
191  if (rank==0)
192  std::cout << "Input adapter for Tpetra::RowGraph" << std::endl;
193 
194  RCP<const ztrowgraph_t> ctrG = rcp_const_cast<const ztrowgraph_t>(
195  rcp_dynamic_cast<ztrowgraph_t>(tG));
196 
197  RCP<adapter_t> trGInput;
198 
199  try {
200  trGInput = rcp(new adapter_t(ctrG));
201  }
202  catch (std::exception &e){
203  aok = false;
204  std::cout << e.what() << std::endl;
205  }
206  TEST_FAIL_AND_EXIT(*comm, aok, "TpetraRowGraphAdapter ", 1);
207 
208  fail = verifyInputAdapter<ztrowgraph_t>(*trGInput, *trG);
209 
210  gfail = globalFail(comm, fail);
211 
212  if (!gfail){
213  ztrowgraph_t *mMigrate = NULL;
214  try{
215  trGInput->applyPartitioningSolution( *trG, mMigrate, solution);
216  newG = rcp(mMigrate);
217  }
218  catch (std::exception &e){
219  fail = 11;
220  }
221 
222  gfail = globalFail(comm, fail);
223 
224  if (!gfail){
225  RCP<const ztrowgraph_t> cnewG =
226  rcp_const_cast<const ztrowgraph_t>(newG);
227  RCP<adapter_t> newInput;
228  try{
229  newInput = rcp(new adapter_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, "TpetraRowGraphAdapter 2 ", 1);
236 
237  if (rank==0){
238  std::cout <<
239  "Input adapter for Tpetra::RowGraph migrated to proc 0" <<
240  std::endl;
241  }
242  fail = verifyInputAdapter<ztrowgraph_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  // DONE
254 
255  if (rank==0)
256  std::cout << "PASS" << std::endl;
257 }
int globalFail(const RCP< const Comm< int > > &comm, int fail)
size_t getLocalNumEdges() const
Returns the number of edges on this process.
Defines TpetraRowGraphAdapter class.
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
int zlno_t
common code used by tests
Tpetra::RowGraph< zlno_t, zgno_t, znode_t > ztrowgraph_t
A PartitioningSolution is a solution to a partitioning problem.
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > ztcrsgraph_t
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::TpetraRowGraphAdapter< User > &ia, ztrowgraph_t &graph)
void printFailureCode(const RCP< const Comm< int > > &comm, int fail)
int main(int argc, char *argv[])
size_t getLocalNumVertices() const
Returns the number of vertices on this process.
void printGraph(RCP< const Comm< int > > &comm, zlno_t nvtx, const zgno_t *vtxIds, const zlno_t *offsets, const zgno_t *edgeIds)
void getVertexIDsView(const gno_t *&ids) const
Sets pointers to this process&#39; graph entries.
void getEdgesView(const lno_t *&offsets, const gno_t *&adjIds) const
Gets adjacency lists for all vertices in a compressed sparse row (CSR) format.
Provides access for Zoltan2 to Tpetra::RowGraph data.
std::string testDataFilePath(".")