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