Zoltan2
Loading...
Searching...
No Matches
XpetraMultiVectorInput.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// @HEADER
46// ***********************************************************************
47// Zoltan2: Sandia Partitioning Ordering & Coloring Library
48// ***********************************************************************
49
55#include <string>
56
60
61#include <Teuchos_DefaultComm.hpp>
62#include <Teuchos_RCP.hpp>
63#include <Teuchos_Comm.hpp>
64#include <Teuchos_CommHelpers.hpp>
65
66using Teuchos::RCP;
67using Teuchos::rcp;
68using Teuchos::rcp_const_cast;
69using Teuchos::Comm;
70
71typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> tvector_t;
72typedef Xpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> xvector_t;
73
75template <typename User>
78{
79 RCP<const Comm<int> > comm = vector.getMap()->getComm();
80 int fail = 0, gfail=0;
81
82 if (!fail && ia.getNumEntriesPerID() !=nvec)
83 fail = 42;
84
85 size_t length = vector.getLocalLength();
86
87 if (!fail && ia.getLocalNumIDs() != length)
88 fail = 4;
89
90 gfail = globalFail(*comm, fail);
91
92 if (!gfail){
93 const zgno_t *vtxIds=NULL;
94 const zscalar_t *vals=NULL;
95 int stride;
96
97 size_t nvals = ia.getLocalNumIDs();
98 if (nvals != vector.getLocalLength())
99 fail = 8;
100
101 ia.getIDsView(vtxIds);
102
103 for (int v=0; v < nvec; v++){
104 auto vecdata = vector.getData(v);
105
106 ia.getEntriesView(vals, stride, v);
107
108 if (!fail && stride != 1)
109 fail = 10;
110
111 // check the values returned
112 for (size_t i = 0; i < length; i++)
113 if (vals[i*stride] != vecdata[i]) {
114 fail = 104;
115 }
116 }
117
118 gfail = globalFail(*comm, fail);
119 }
120
121
122 return gfail;
123}
124
125
126template <typename User>
129 std::vector<const zscalar_t *> &weights, std::vector<int> &strides)
130{
131 // Check the input adapter
132 int gfail = verifyInputAdapter(ia, vector, nvec);
133 int fail = gfail;
134
135 RCP<const Comm<int> > comm = vector.getMap()->getComm();
136 int wdim = weights.size();
137
138 // Check the input adapter weights
139 if (!gfail && (ia.getNumWeightsPerID() != wdim)) fail = 41;
140
141 gfail = globalFail(*comm, fail);
142
143 if (!gfail && wdim) {
144 const zscalar_t *wgt =NULL;
145 int stride;
146
147 for (int w=0; !fail && w < wdim; w++){
148 ia.getWeightsView(wgt, stride, w);
149
150 if (!fail && stride != strides[w])
151 fail = 101;
152
153 for (size_t v=0; !fail && v < vector.getLocalLength(); v++){
154 if (wgt[v*stride] != weights[w][v*stride])
155 fail=102;
156 }
157 }
158
159 gfail = globalFail(*comm, fail);
160 }
161 return gfail;
162}
163
164
166
168 std::ifstream &fp,
169 size_t &nIDs,
170 size_t &nEdges,
171 char *code,
172 int &nWgts)
173{
174 // Read the header info from a Chaco .graph file
175 std::string line;
176 std::getline(fp, line);
177 while (line[0]=='#') std::getline(fp, line); // skip comments
178 std::istringstream issHeader(line);
179 issHeader >> nIDs >> nEdges >> code;
180 if (!strcmp(code, "010") || !strcmp(code, "011")) {
181 if (!(issHeader >> nWgts)) nWgts = 1;
182 }
183}
184
185template <typename User>
188 const char *fileprefixInp,
189 const Teuchos::Comm<int> &comm
190)
191{
192 // Compares files generated by generateFiles (in fileprefixGen.*)
193 // to input files (in fileprefixInp.*).
194
195 int fail = 0, gfail=0;
196
197 std::string tmp(fileprefixInp);
198 tmp = tmp + "_generated";
199 const char *fileprefixGen = tmp.c_str();
200
201 ia.generateFiles(fileprefixGen, comm);
202
203 // Only rank zero needs to check the resulting files
204 if (comm.getRank() == 0) {
205
206 size_t nIDsGen, nIDsInp;
207 size_t nEdgesGen, nEdgesInp;
208 char codeGen[4], codeInp[4];
209 int nWgtsGen = 0, nWgtsInp = 0;
210 std::string lineGen, lineInp;
211
212 std::ifstream fpGen, fpInp;
213 std::string graphFilenameGen = fileprefixGen;
214 graphFilenameGen = graphFilenameGen + ".graph";
215 std::string graphFilenameInp = fileprefixInp;
216 graphFilenameInp = graphFilenameInp + ".graph";
217
218 // Read header info from generated file
219 fpGen.open(graphFilenameGen.c_str(), std::ios::in);
220 readChacoGraphHeaderInfo(fpGen, nIDsGen, nEdgesGen, codeGen, nWgtsGen);
221
222 // Read header info from input file
223 fpInp.open(graphFilenameInp.c_str(), std::ios::in);
224 readChacoGraphHeaderInfo(fpInp, nIDsInp, nEdgesInp, codeInp, nWgtsInp);
225
226 // input file and generated file should have same number of IDs
227 if (nIDsGen != nIDsInp) {
228 std::cout << "GenerateFiles: nIDsGen " << nIDsGen
229 << " != nIDsInp " << nIDsInp << std::endl;
230 fail = 2222;
231 }
232
233 // Vector adapters don't have edges
234 if (!fail && nEdgesGen != 0) {
235 std::cout << "GenerateFiles: nEdgesGen " << nEdgesGen << " != 0"
236 << std::endl;
237 fail = 2223;
238 }
239
240 // Check the weights, if any
241 if (!fail && nWgtsGen) {
242 if (nWgtsInp) {
243 // input file has weights; compare weights
244 size_t cntWgtLines;
245 for (cntWgtLines = 0; cntWgtLines < nIDsGen &&
246 std::getline(fpGen, lineGen) &&
247 std::getline(fpInp, lineInp); cntWgtLines++) {
248
249 std::istringstream issGen(lineGen);
250 std::istringstream issInp(lineInp);
251
252 int nw = 0;
253 double wgtGen, wgtInp;
254
255 while (issGen >> wgtGen) {
256
257 if (nw < nWgtsInp) {
258 issInp >> wgtInp;
259 if (wgtGen != wgtInp) fail = 2224; // weights don't match
260 }
261 nw++;
262 }
263
264 if (nw != nWgtsGen) fail = 2225; // Too many or few weights on line
265 }
266 if (cntWgtLines != nIDsGen) fail = 2226; // not enough input lines
267 }
268 else {
269 // input file does not have weights;
270 // just make sure generated file has enough weights
271 size_t cntWgtLines = 0;
272 for (cntWgtLines = 0; cntWgtLines < nIDsGen &&
273 std::getline(fpGen, lineGen); cntWgtLines++) {
274 std::istringstream issGen(lineGen);
275 int nw = 0;
276 double wgtGen;
277 while (issGen) { issGen >> wgtGen; nw++; }
278 if (nw != nWgtsGen) fail = 2227; // Too many or few weights on line
279 }
280 if (cntWgtLines != nIDsGen) fail = 2228; // not enough input lines
281 }
282 }
283
284 fpGen.close();
285 fpInp.close();
286
287 // check coordinate files
288 if (!fail) {
289 std::string coordsFilenameGen = fileprefixGen;
290 coordsFilenameGen = coordsFilenameGen + ".coords";
291 std::string coordsFilenameInp = fileprefixInp;
292 coordsFilenameInp = coordsFilenameInp + ".coords";
293
294 fpGen.open(coordsFilenameGen.c_str(), std::ios::in);
295 fpInp.open(coordsFilenameInp.c_str(), std::ios::in);
296
297 size_t cnt;
298 for (cnt = 0; std::getline(fpGen,lineGen) &&
299 std::getline(fpInp,lineInp); cnt++) {
300
301 // Check each token
302 std::istringstream issGen(lineGen);
303 std::istringstream issInp(lineInp);
304
305 double xGen, xInp;
306 issGen >> xGen;
307 issInp >> xInp;
308 while (issGen && issInp) {
309
310 if (xGen != xInp) {
311 std::cout << "Coordinates " << xGen << " != " << xInp
312 << std::endl;
313 fail = 333;
314 }
315 issGen >> xGen;
316 issInp >> xInp;
317 }
318
319 // Check same number of tokens: are there any left in either line?
320 if (issGen || issInp) {
321 std::cout << "Dimension of generated file != dimension of input file"
322 << std::endl;
323 fail = 334;
324 }
325 }
326
327 // Did we have the correct number of coordinates?
328 if (!fail && cnt != nIDsGen) {
329 std::cout << "Number of coordinates read " << cnt
330 << " != number of IDs " << nIDsGen << std::endl;
331 fail = 444;
332 }
333
334 fpGen.close();
335 fpInp.close();
336 }
337 }
338
339 gfail = globalFail(comm, fail);
340 return gfail;
341}
342
344
345int main(int narg, char *arg[])
346{
347 Tpetra::ScopeGuard tscope(&narg, &arg);
348 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
349
350 int rank = comm->getRank();
351 int fail = 0, gfail=0;
352 bool aok = true;
353
354 // Read run-time options.
355 std::string inputFilePrefix("simple");
356 Teuchos::CommandLineProcessor cmdp (false, false);
357 cmdp.setOption("file", &inputFilePrefix,
358 "chaco file (prefix) to be used in test");
359 cmdp.parse(narg, arg);
360
361 // Create object that can give us test Tpetra, Xpetra
362 // and Epetra vectors for testing.
363
364 RCP<UserInputForTests> uinput;
365 Teuchos::ParameterList params;
366 params.set("input file", inputFilePrefix);
367 params.set("file type", "Chaco");
368
369 try{
370 uinput = rcp(new UserInputForTests(params, comm));
371 }
372 catch(std::exception &e){
373 aok = false;
374 std::cout << e.what() << std::endl;
375 }
376 TEST_FAIL_AND_EXIT(*comm, aok, "input ", 1);
377
378 RCP<tvector_t> inputMVector; // original vector (for checking)
379 RCP<tvector_t> migratedMVector; // migrated vector
380
381 int nVec = 2;
382
383 inputMVector = uinput->getUICoordinates();
384 size_t vlen = inputMVector->getLocalLength();
385
386 // Vertex weights
387 std::vector<const zscalar_t *> IDWeights;
388 std::vector<int> IDWeightsStrides;
389
390 int nWeights = 0;
391 if (uinput->hasUIWeights()) {
392
393 auto uiweights = uinput->getUIWeights();
394
395 nWeights = uiweights->getNumVectors();
396 IDWeights.resize(nWeights);
397 IDWeightsStrides.resize(nWeights);
398
399 for (int w = 0; w < nWeights; w++) {
400 IDWeights[w] = uiweights->getData(w).getRawPtr();
401 IDWeightsStrides[w] = 1;
402 }
403 }
404
405 // To test migration in the input adapter we need a Solution object.
406
407 RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment(comm));
408
411 typedef ia_t::part_t part_t;
412
413 part_t *p = new part_t [vlen];
414 memset(p, 0, sizeof(part_t) * vlen);
415 ArrayRCP<part_t> solnParts(p, 0, vlen, true);
416
417 soln_t solution(env, comm, nWeights);
418 solution.setParts(solnParts);
419
420
422 // User object is Tpetra::MultiVector, no weights
423 if (!gfail){
424 if (rank==0)
425 std::cout << "Constructed with Tpetra::MultiVector" << std::endl;
426
427 RCP<const tvector_t> ctV = rcp_const_cast<const tvector_t>(inputMVector);
428 RCP<Zoltan2::XpetraMultiVectorAdapter<tvector_t> > tVInput;
429
430 try {
431 tVInput =
433 IDWeights, IDWeightsStrides));
434 }
435 catch (std::exception &e){
436 aok = false;
437 std::cout << e.what() << std::endl;
438 }
439 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter ", 1);
440
441 fail = verifyInputAdapter<tvector_t>(*tVInput, *inputMVector, nVec,
442 IDWeights, IDWeightsStrides);
443 fail = verifyGenerateFiles(*tVInput, inputFilePrefix.c_str(), *comm);
444
445 gfail = globalFail(*comm, fail);
446
447 if (!gfail){
448 tvector_t *vMigrate = NULL;
449 try{
450 tVInput->applyPartitioningSolution(*inputMVector, vMigrate, solution);
451 migratedMVector = rcp(vMigrate);
452 }
453 catch (std::exception &e){
454 fail = 11;
455 }
456
457 gfail = globalFail(*comm, fail);
458
459 if (!gfail){
460 RCP<const tvector_t> cnewV =
461 rcp_const_cast<const tvector_t>(migratedMVector);
462 RCP<Zoltan2::XpetraMultiVectorAdapter<tvector_t> > newInput;
463 try{
464 newInput =
466 }
467 catch (std::exception &e){
468 aok = false;
469 std::cout << e.what() << std::endl;
470 }
471 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 2 ", 1);
472
473 if (rank==0){
474 std::cout << "Constructed with ";
475 std::cout << "Tpetra::MultiVector migrated to proc 0" << std::endl;
476 }
477 // We didn't apply partitioning soln to weights,
478 // so don't check them when verifying adapter
479 fail = verifyInputAdapter<tvector_t>(*newInput, *migratedMVector, nVec);
480 if (fail) fail += 100;
481 gfail = globalFail(*comm, fail);
482 }
483 }
484 if (gfail){
485 printFailureCode(*comm, fail);
486 }
487 }
488
490 // User object is Xpetra::MultiVector
491 if (!gfail){
492 if (rank==0)
493 std::cout << "Constructed with Xpetra::MultiVector" << std::endl;
494
495 RCP<xvector_t> xV =
497 RCP<const xvector_t> cxV = rcp_const_cast<const xvector_t>(xV);
498 RCP<Zoltan2::XpetraMultiVectorAdapter<xvector_t> > xVInput;
499
500 try {
501 xVInput =
503 IDWeights,
504 IDWeightsStrides));
505 }
506 catch (std::exception &e){
507 aok = false;
508 std::cout << e.what() << std::endl;
509 }
510 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 3 ", 1);
511
512 fail = verifyInputAdapter<xvector_t>(*xVInput, *inputMVector, nVec,
513 IDWeights, IDWeightsStrides);
514
515 gfail = globalFail(*comm, fail);
516
517 if (!gfail){
518 xvector_t *vMigrate =NULL;
519 try{
520 xVInput->applyPartitioningSolution(*xV, vMigrate, solution);
521 }
522 catch (std::exception &e){
523 fail = 11;
524 }
525
526 gfail = globalFail(*comm, fail);
527
528 if (!gfail){
529 RCP<const xvector_t> cnewV(vMigrate);
530 RCP<Zoltan2::XpetraMultiVectorAdapter<xvector_t> > newInput;
531 try{
532 newInput =
534 }
535 catch (std::exception &e){
536 aok = false;
537 std::cout << e.what() << std::endl;
538 }
539 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 4 ", 1);
540
541 if (rank==0){
542 std::cout << "Constructed with ";
543 std::cout << "Xpetra::MultiVector migrated to proc 0" << std::endl;
544 }
545 fail = verifyInputAdapter<xvector_t>(*newInput, *migratedMVector, nVec);
546 if (fail) fail += 100;
547 gfail = globalFail(*comm, fail);
548 }
549 }
550 if (gfail){
551 printFailureCode(*comm, fail);
552 }
553 }
554
555#ifdef HAVE_EPETRA_DATA_TYPES
557 // User object is Epetra_MultiVector
558 typedef Epetra_MultiVector evector_t;
559 if (!gfail){
560 if (rank==0)
561 std::cout << "Constructed with Epetra_MultiVector" << std::endl;
562
563 RCP<evector_t> eV =
564 rcp(new Epetra_MultiVector(uinput->getUIEpetraCrsGraph()->RowMap(),
565 nVec));
566 for (int v = 0; v < nVec; v++) {
567 auto inV = inputMVector->getData(v);
568 for (int i = 0; i < eV->MyLength(); i++)
569 eV->ReplaceMyValue(i, v, inV[i]);
570 }
571
572 RCP<const evector_t> ceV = rcp_const_cast<const evector_t>(eV);
573 RCP<Zoltan2::XpetraMultiVectorAdapter<evector_t> > eVInput;
574
575 bool goodAdapter = true;
576 try {
577 eVInput =
579 IDWeights, IDWeightsStrides));
580 }
581 catch (std::exception &e){
582 if (std::is_same<znode_t, Xpetra::EpetraNode>::value) {
583 aok = false;
584 goodAdapter = false;
585 std::cout << e.what() << std::endl;
586 }
587 else {
588 // We expect an error from Xpetra when znode_t != Xpetra::EpetraNode
589 // Ignore it, but skip tests using this matrix adapter.
590 std::cout << "Node type is not supported by Xpetra's Epetra interface;"
591 << " Skipping this test." << std::endl;
592 std::cout << "FYI: Here's the exception message: " << std::endl
593 << e.what() << std::endl;
594 goodAdapter = false;
595 }
596 }
597 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 5 ", 1);
598
599 if (goodAdapter) {
600 fail = verifyInputAdapter<evector_t>(*eVInput, *inputMVector, nVec,
601 IDWeights, IDWeightsStrides);
602
603 gfail = globalFail(*comm, fail);
604
605 if (!gfail){
606 evector_t *vMigrate =NULL;
607 try{
608 eVInput->applyPartitioningSolution(*eV, vMigrate, solution);
609 }
610 catch (std::exception &e){
611 fail = 11;
612 }
613
614 gfail = globalFail(*comm, fail);
615
616 if (!gfail){
617 RCP<const evector_t> cnewV(vMigrate, true);
618 RCP<Zoltan2::XpetraMultiVectorAdapter<evector_t> > newInput;
619 try{
620 newInput =
622 }
623 catch (std::exception &e){
624 aok = false;
625 std::cout << e.what() << std::endl;
626 }
627 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 6 ", 1);
628
629 if (rank==0){
630 std::cout << "Constructed with ";
631 std::cout << "Epetra_MultiVector migrated to proc 0" << std::endl;
632 }
633 fail = verifyInputAdapter<evector_t>(*newInput, *migratedMVector,
634 nVec);
635 if (fail) fail += 100;
636 gfail = globalFail(*comm, fail);
637 }
638 }
639 if (gfail){
640 printFailureCode(*comm, fail);
641 }
642 }
643 }
644#endif
645
647 // DONE
648
649 if (rank==0)
650 std::cout << "PASS" << std::endl;
651}
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)
int verifyInputAdapter(Zoltan2::XpetraMultiVectorAdapter< User > &ia, tvector_t &vector, int nvec)
void readChacoGraphHeaderInfo(std::ifstream &fp, size_t &nIDs, size_t &nEdges, char *code, int &nWgts)
Xpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > xvector_t
int verifyGenerateFiles(Zoltan2::VectorAdapter< User > &ia, const char *fileprefixInp, const Teuchos::Comm< int > &comm)
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tvector_t
Traits for application input objects.
common code used by tests
float zscalar_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.
VectorAdapter defines the interface for vector input.
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.