Zoltan2
Loading...
Searching...
No Matches
XpetraVectorInput.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
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;
67
68typedef Tpetra::Vector<zscalar_t, zlno_t, zgno_t, znode_t> tvector_t;
69typedef Xpetra::Vector<zscalar_t, zlno_t, zgno_t, znode_t> xvector_t;
70
71void printVector(RCP<const Comm<int> > &comm, zlno_t vlen,
72 const zgno_t *vtxIds, const zscalar_t *vals)
73{
74 int rank = comm->getRank();
75 int nprocs = comm->getSize();
76 comm->barrier();
77 for (int p=0; p < nprocs; p++){
78 if (p == rank){
79 std::cout << rank << ":" << std::endl;
80 for (zlno_t i=0; i < vlen; i++){
81 std::cout << " " << vtxIds[i] << ": " << vals[i] << std::endl;
82 }
83 std::cout.flush();
84 }
85 comm->barrier();
86 }
87 comm->barrier();
88}
89
90template <typename User>
93 zscalar_t **weights, int *strides)
94{
95 RCP<const Comm<int> > comm = vector.getMap()->getComm();
96 int fail = 0, gfail=0;
97
98 if (!fail && ia.getNumEntriesPerID() !=1)
99 fail = 42;
100
101 if (!fail && ia.getNumWeightsPerID() !=wdim)
102 fail = 41;
103
104 if (!fail && ia.getLocalNumIDs() != vector.getLocalLength())
105 fail = 4;
106
107 gfail = globalFail(*comm, fail);
108
109 if (!gfail){
110 const zgno_t *vtxIds=NULL;
111 const zscalar_t *vals=NULL;
112 int stride;
113
114 size_t nvals = ia.getLocalNumIDs();
115 if (nvals != vector.getLocalLength())
116 fail = 8;
117
118 ia.getIDsView(vtxIds);
119 ia.getEntriesView(vals, stride);
120
121 if (!fail && stride != 1)
122 fail = 10;
123
124 gfail = globalFail(*comm, fail);
125
126 if (gfail == 0){
127 printVector(comm, nvals, vtxIds, vals);
128 }
129 }
130
131 if (!gfail && wdim){
132 const zscalar_t *wgt =NULL;
133 int stride;
134
135 for (int w=0; !fail && w < wdim; w++){
136 ia.getWeightsView(wgt, stride, w);
137
138 if (!fail && stride != strides[w])
139 fail = 101;
140
141 for (size_t v=0; !fail && v < vector.getLocalLength(); v++){
142 if (wgt[v*stride] != weights[w][v*stride])
143 fail=102;
144 }
145 }
146
147 gfail = globalFail(*comm, fail);
148 }
149
150 return gfail;
151}
152
153
154int main(int narg, char *arg[])
155{
156 Tpetra::ScopeGuard tscope(&narg, &arg);
157 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
158
159 int rank = comm->getRank();
160 int fail = 0, gfail=0;
161 bool aok = true;
162
163 // Create object that can give us test Tpetra, Xpetra
164 // and Epetra vectors for testing.
165
166 RCP<UserInputForTests> uinput;
167 Teuchos::ParameterList params;
168 params.set("input file", "simple");
169 params.set("file type", "Chaco");
170 try{
171 uinput = rcp(new UserInputForTests(params, comm));
172 }
173 catch(std::exception &e){
174 aok = false;
175 std::cout << e.what() << std::endl;
176 }
177 TEST_FAIL_AND_EXIT(*comm, aok, "input ", 1);
178
179 RCP<tvector_t> tV; // original vector (for checking)
180 RCP<tvector_t> newV; // migrated vector
181
182 tV = rcp(new tvector_t(uinput->getUITpetraCrsGraph()->getRowMap()));
183 tV->randomize();
184 size_t vlen = tV->getLocalLength();
185
186 // To test migration in the input adapter we need a Solution object.
187
188 RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment(comm));
189
190 int nWeights = 1;
191
193 typedef adapter_t::part_t part_t;
194
195 part_t *p = new part_t [vlen];
196 memset(p, 0, sizeof(part_t) * vlen);
197 ArrayRCP<part_t> solnParts(p, 0, vlen, true);
198
199 std::vector<const zscalar_t *> emptyWeights;
200 std::vector<int> emptyStrides;
201
202 Zoltan2::PartitioningSolution<adapter_t> solution(env, comm, nWeights);
203 solution.setParts(solnParts);
204
206 // User object is Tpetra::Vector, no weights
207 if (!gfail){
208 if (rank==0)
209 std::cout << "Constructed with Tpetra::Vector" << std::endl;
210
211 RCP<const tvector_t> ctV = rcp_const_cast<const tvector_t>(tV);
212 RCP<adapter_t> tVInput;
213
214 try {
215 tVInput =
217 emptyWeights, emptyStrides));
218 }
219 catch (std::exception &e){
220 aok = false;
221 std::cout << e.what() << std::endl;
222 }
223 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter ", 1);
224
225 fail = verifyInputAdapter<tvector_t>(*tVInput, *tV, 0, NULL, NULL);
226
227 gfail = globalFail(*comm, fail);
228
229 if (!gfail){
230 tvector_t *vMigrate = NULL;
231 try{
232 tVInput->applyPartitioningSolution(*tV, vMigrate, solution);
233 newV = rcp(vMigrate);
234 }
235 catch (std::exception &e){
236 fail = 11;
237 }
238
239 gfail = globalFail(*comm, fail);
240
241 if (!gfail){
242 RCP<const tvector_t> cnewV = rcp_const_cast<const tvector_t>(newV);
243 RCP<Zoltan2::XpetraMultiVectorAdapter<tvector_t> > newInput;
244 try{
245 newInput = rcp(new Zoltan2::XpetraMultiVectorAdapter<tvector_t>(cnewV,
246 emptyWeights, emptyStrides));
247 }
248 catch (std::exception &e){
249 aok = false;
250 std::cout << e.what() << std::endl;
251 }
252 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 2 ", 1);
253
254 if (rank==0){
255 std::cout << "Constructed with ";
256 std::cout << "Tpetra::Vector migrated to proc 0" << std::endl;
257 }
258 fail = verifyInputAdapter<tvector_t>(*newInput, *newV, 0, NULL, NULL);
259 if (fail) fail += 100;
260 gfail = globalFail(*comm, fail);
261 }
262 }
263 if (gfail){
264 printFailureCode(*comm, fail);
265 }
266 }
267
269 // User object is Xpetra::Vector
270 if (!gfail){
271 if (rank==0)
272 std::cout << "Constructed with Xpetra::Vector" << std::endl;
273
274 RCP<tvector_t> xtV =
275 rcp(new tvector_t(uinput->getUITpetraCrsGraph()->getRowMap()));
276 xtV->randomize();
278 RCP<const xvector_t> cxV = rcp_const_cast<const xvector_t>(xV);
279 RCP<Zoltan2::XpetraMultiVectorAdapter<xvector_t> > xVInput;
280
281 try {
282 xVInput =
284 emptyWeights, emptyStrides));
285 }
286 catch (std::exception &e){
287 aok = false;
288 std::cout << e.what() << std::endl;
289 }
290 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 3 ", 1);
291
292 fail = verifyInputAdapter<xvector_t>(*xVInput, *tV, 0, NULL, NULL);
293
294 gfail = globalFail(*comm, fail);
295
296 if (!gfail){
297 xvector_t *vMigrate =NULL;
298 try{
299 xVInput->applyPartitioningSolution(*xV, vMigrate, solution);
300 }
301 catch (std::exception &e){
302 fail = 11;
303 }
304
305 gfail = globalFail(*comm, fail);
306
307 if (!gfail){
308 RCP<const xvector_t> cnewV(vMigrate);
309 RCP<Zoltan2::XpetraMultiVectorAdapter<xvector_t> > newInput;
310 try{
311 newInput =
313 emptyWeights, emptyStrides));
314 }
315 catch (std::exception &e){
316 aok = false;
317 std::cout << e.what() << std::endl;
318 }
319 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 4 ", 1);
320
321 if (rank==0){
322 std::cout << "Constructed with ";
323 std::cout << "Xpetra::Vector migrated to proc 0" << std::endl;
324 }
325 fail = verifyInputAdapter<xvector_t>(*newInput, *newV, 0, NULL, NULL);
326 if (fail) fail += 100;
327 gfail = globalFail(*comm, fail);
328 }
329 }
330 if (gfail){
331 printFailureCode(*comm, fail);
332 }
333 }
334
335#ifdef HAVE_EPETRA_DATA_TYPES
337 // User object is Epetra_Vector
338 typedef Epetra_Vector evector_t;
339 if (!gfail){
340 if (rank==0)
341 std::cout << "Constructed with Epetra_Vector" << std::endl;
342
343 RCP<evector_t> eV =
344 rcp(new Epetra_Vector(uinput->getUIEpetraCrsGraph()->RowMap()));
345 eV->Random();
346 RCP<const evector_t> ceV = rcp_const_cast<const evector_t>(eV);
347 RCP<Zoltan2::XpetraMultiVectorAdapter<evector_t> > eVInput;
348
349 bool goodAdapter = true;
350 try {
351 eVInput =
353 emptyWeights, emptyStrides));
354 }
355 catch (std::exception &e){
356 if (std::is_same<znode_t, Xpetra::EpetraNode>::value) {
357 aok = false;
358 goodAdapter = false;
359 std::cout << e.what() << std::endl;
360 }
361 else {
362 // We expect an error from Xpetra when znode_t != Xpetra::EpetraNode
363 // Ignore it, but skip tests using this matrix adapter.
364 std::cout << "Node type is not supported by Xpetra's Epetra interface;"
365 << " Skipping this test." << std::endl;
366 std::cout << "FYI: Here's the exception message: " << std::endl
367 << e.what() << std::endl;
368 goodAdapter = false;
369 }
370 }
371 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 5 ", 1);
372
373 if (goodAdapter) {
374 fail = verifyInputAdapter<evector_t>(*eVInput, *tV, 0, NULL, NULL);
375
376 gfail = globalFail(*comm, fail);
377
378 if (!gfail){
379 evector_t *vMigrate =NULL;
380 try{
381 eVInput->applyPartitioningSolution(*eV, vMigrate, solution);
382 }
383 catch (std::exception &e){
384 fail = 11;
385 }
386
387 gfail = globalFail(*comm, fail);
388
389 if (!gfail){
390 RCP<const evector_t> cnewV(vMigrate, true);
391 RCP<Zoltan2::XpetraMultiVectorAdapter<evector_t> > newInput;
392 try{
393 newInput =
395 emptyWeights, emptyStrides));
396 }
397 catch (std::exception &e){
398 aok = false;
399 std::cout << e.what() << std::endl;
400 }
401 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 6 ", 1);
402
403 if (rank==0){
404 std::cout << "Constructed with ";
405 std::cout << "Epetra_Vector migrated to proc 0" << std::endl;
406 }
407 fail = verifyInputAdapter<evector_t>(*newInput, *newV, 0, NULL, NULL);
408 if (fail) fail += 100;
409 gfail = globalFail(*comm, fail);
410 }
411 }
412 if (gfail){
413 printFailureCode(*comm, fail);
414 }
415 }
416 }
417#endif
418
420 // DONE
421
422 if (rank==0)
423 std::cout << "PASS" << std::endl;
424}
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)
Xpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > xvector_t
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tvector_t
Xpetra::Vector< zscalar_t, zlno_t, zgno_t, znode_t > xvector_t
Tpetra::Vector< zscalar_t, zlno_t, zgno_t, znode_t > tvector_t
int verifyInputAdapter(Zoltan2::XpetraMultiVectorAdapter< User > &ia, tvector_t &vector, int wdim, zscalar_t **weights, int *strides)
void printVector(RCP< const Comm< int > > &comm, zlno_t vlen, const zgno_t *vtxIds, const zscalar_t *vals)
Traits for application input objects.
common code used by tests
float zscalar_t
Tpetra::Map ::local_ordinal_type zlno_t
Tpetra::Map ::global_ordinal_type zgno_t
Defines the XpetraMultiVectorAdapter.
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.
void setParts(ArrayRCP< part_t > &partList)
The algorithm uses setParts to set the solution.
void getEntriesView(const scalar_t *&elements, int &stride, int idx=0) const
Provide a pointer to the elements of the specified vector.
int getNumEntriesPerID() const
Return the number of vectors.
int getNumWeightsPerID() const
Returns the number of weights per object. Number of weights per object should be zero or greater....
void getWeightsView(const scalar_t *&weights, int &stride, int idx) const
size_t getLocalNumIDs() const
Returns the number of objects on this process.
static const std::string fail
static ArrayRCP< ArrayRCP< zscalar_t > > weights
SparseMatrixAdapter_t::part_t part_t
static RCP< User > convertToXpetra(const RCP< User > &a)
Convert the object to its Xpetra wrapped version.