Zoltan2
AdapterForTests.hpp
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 
51 #ifndef ADAPTERFORTESTS
52 #define ADAPTERFORTESTS
53 
54 #include <Zoltan2_Parameters.hpp>
55 #include <UserInputForTests.hpp>
56 
59 
65 
66 #ifdef HAVE_ZOLTAN2_PAMGEN
68 #endif
69 
70 #include <Teuchos_DefaultComm.hpp>
71 #include <Teuchos_XMLObject.hpp>
72 #include <Teuchos_FileInputSource.hpp>
73 
74 #include <Tpetra_MultiVector.hpp>
75 #include <Tpetra_CrsMatrix.hpp>
76 
77 #include <string>
78 #include <iostream>
79 #include <vector>
80 
81 using namespace std;
82 using Teuchos::RCP;
83 using Teuchos::ArrayRCP;
84 using Teuchos::ArrayView;
85 using Teuchos::Array;
86 using Teuchos::Comm;
87 using Teuchos::rcp;
88 using Teuchos::arcp;
89 using Teuchos::rcp_const_cast;
90 using Teuchos::ParameterList;
91 using std::string;
92 
93 /* \brief A class for constructing Zoltan2 input adapters */
95 public:
96 
101 
106 
114 
115 #ifdef HAVE_ZOLTAN2_PAMGEN
117 #else
118  // This typedef exists only to satisfy the compiler.
119  // PamgenMeshAdapter cannot be used when Trilinos is not built with Pamgen
121 #endif
122 
132  static base_adapter_t* getAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
133 
134 private:
143  static base_adapter_t*
144  getBasicIdentiferAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
145 
154  static base_adapter_t*
155  getXpetraMVAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
156 
165  static base_adapter_t*
166  getXpetraCrsGraphAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
167 
176  static base_adapter_t*
177  getXpetraCrsMatrixAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
178 
187  static base_adapter_t*
188  getBasicVectorAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
189 
198  static base_adapter_t*
199  getPamgenMeshAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
200 
201 
210  template <typename T>
211  static void InitializeVectorData(const RCP<T> &data,
212  vector<const zscalar_t *> &coords,
213  vector<int> & strides,
214  int stride);
215 
216 #ifdef HAVE_EPETRA_DATA_TYPES
217 
225  template <typename T>
226  static void InitializeEpetraVectorData(const RCP<T> &data,
227  vector<const zscalar_t *> &coords,
228  vector<int> & strides,
229  int stride);
230 #endif
231 };
232 
233 
235  const ParameterList &pList,
236  const RCP<const Comm<int> > &comm)
237 {
238  AdapterForTests::base_adapter_t * ia = nullptr; // input adapter
239 
240  if(!pList.isParameter("input adapter"))
241  {
242  std::cerr << "Input adapter unspecified" << std::endl;
243  return ia;
244  }
245 
246  // pick method for chosen adapter
247  string adapter_name = pList.get<string>("input adapter");
248  if(adapter_name == "BasicIdentifier")
249  ia = AdapterForTests::getBasicIdentiferAdapterForInput(uinput, pList, comm);
250  else if(adapter_name == "XpetraMultiVector")
251  ia = AdapterForTests::getXpetraMVAdapterForInput(uinput, pList, comm);
252  else if(adapter_name == "XpetraCrsGraph")
253  ia = getXpetraCrsGraphAdapterForInput(uinput,pList, comm);
254  else if(adapter_name == "XpetraCrsMatrix")
255  ia = getXpetraCrsMatrixAdapterForInput(uinput,pList, comm);
256  else if(adapter_name == "BasicVector")
257  ia = getBasicVectorAdapterForInput(uinput,pList, comm);
258  else if(adapter_name == "PamgenMesh")
259  ia = getPamgenMeshAdapterForInput(uinput,pList, comm);
260  else
261  std::cerr << "Input adapter type: " + adapter_name + ", is unavailable, or misspelled." << std::endl;
262 
263  return ia;
264 }
265 
266 
267 AdapterForTests::base_adapter_t * AdapterForTests::getBasicIdentiferAdapterForInput(UserInputForTests *uinput,
268  const ParameterList &pList,
269  const RCP<const Comm<int> > &comm)
270 {
271 
272  if(!pList.isParameter("data type"))
273  {
274  std::cerr << "Input data type unspecified" << std::endl;
275  return nullptr;
276  }
277 
278  string input_type = pList.get<string>("data type"); // get the input type
279 
280  if (!uinput->hasInputDataType(input_type))
281  {
282  std::cerr << "Input type: " + input_type + " unavailable or misspelled."
283  << std::endl; // bad type
284  return nullptr;
285  }
286 
287  vector<const zscalar_t *> weights;
288  std::vector<int> weightStrides;
289  const zgno_t * globalIds;
290  size_t localCount = 0;
291 
292  // get weights if any
293  // get weights if any
294  if(uinput->hasUIWeights())
295  {
296  RCP<tMVector_t> vtx_weights = uinput->getUIWeights();
297  // copy to weight
298  size_t cols = vtx_weights->getNumVectors();
299  for (size_t i = 0; i< cols; i++) {
300  weights.push_back(vtx_weights->getData(i).getRawPtr());
301  weightStrides.push_back((int)vtx_weights->getStride());
302  }
303  }
304 
305  if(input_type == "coordinates")
306  {
307  RCP<tMVector_t> data = uinput->getUICoordinates();
308  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
309  localCount = data->getLocalLength();
310  }
311  else if(input_type == "tpetra_vector")
312  {
313  RCP<tVector_t> data = uinput->getUITpetraVector();
314  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
315  localCount = data->getLocalLength();
316  }
317  else if(input_type == "tpetra_multivector")
318  {
319  int nvec = pList.get<int>("vector_dimension");
320  RCP<tMVector_t> data = uinput->getUITpetraMultiVector(nvec);
321  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
322  localCount = data->getLocalLength();
323  }
324  else if(input_type == "tpetra_crs_graph")
325  {
326  RCP<tcrsGraph_t> data = uinput->getUITpetraCrsGraph();
327  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
328  localCount = data->getNodeNumCols();
329  }
330  else if(input_type == "tpetra_crs_matrix")
331  {
332  RCP<tcrsMatrix_t> data = uinput->getUITpetraCrsMatrix();
333  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
334  localCount = data->getNodeNumCols();
335  }
336  else if(input_type == "xpetra_vector")
337  {
338  RCP<xVector_t> data = uinput->getUIXpetraVector();
339  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
340  localCount = data->getLocalLength();
341  }
342  else if(input_type == "xpetra_multivector")
343  {
344  int nvec = pList.get<int>("vector_dimension");
345  RCP<xMVector_t> data = uinput->getUIXpetraMultiVector(nvec);
346  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
347  localCount = data->getLocalLength();
348  }
349  else if(input_type == "xpetra_crs_graph")
350  {
351  RCP<xcrsGraph_t> data = uinput->getUIXpetraCrsGraph();
352  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
353  localCount = data->getNodeNumCols();
354  }
355  else if(input_type == "xpetra_crs_matrix")
356  {
357  RCP<xcrsMatrix_t> data = uinput->getUIXpetraCrsMatrix();
358  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
359  localCount = data->getNodeNumCols();
360  }
361 #ifdef HAVE_EPETRA_DATA_TYPES
362  else if(input_type == "epetra_vector")
363  {
364  RCP<Epetra_Vector> data = uinput->getUIEpetraVector();
365  globalIds = (zgno_t *)data->Map().MyGlobalElements();
366  localCount = data->MyLength();
367  }
368  else if(input_type == "epetra_multivector")
369  {
370  int nvec = pList.get<int>("vector_dimension");
371  RCP<Epetra_MultiVector> data = uinput->getUIEpetraMultiVector(nvec);
372  globalIds = (zgno_t *)data->Map().MyGlobalElements();
373  localCount = data->MyLength();
374  }
375  else if(input_type == "epetra_crs_graph")
376  {
377  RCP<Epetra_CrsGraph> data = uinput->getUIEpetraCrsGraph();
378  globalIds = (zgno_t *)data->Map().MyGlobalElements();
379  localCount = data->NumMyCols();
380  }
381  else if(input_type == "epetra_crs_matrix")
382  {
383  RCP<Epetra_CrsMatrix> data = uinput->getUIEpetraCrsMatrix();
384  globalIds = (zgno_t *)data->Map().MyGlobalElements();
385  localCount = data->NumMyCols();
386  }
387 #endif
388 
389  if(localCount == 0) return nullptr;
390  return reinterpret_cast<AdapterForTests::base_adapter_t *>( new AdapterForTests::basic_id_t(zlno_t(localCount),globalIds,weights,weightStrides));
391 }
392 
393 
394 AdapterForTests::base_adapter_t * AdapterForTests::getXpetraMVAdapterForInput(
395  UserInputForTests *uinput,
396  const ParameterList &pList,
397  const RCP<const Comm<int> > &comm)
398 {
399  AdapterForTests::base_adapter_t * adapter = nullptr;
400 
401  if(!pList.isParameter("data type"))
402  {
403  std::cerr << "Input data type unspecified" << std::endl;
404  return adapter;
405  }
406 
407  string input_type = pList.get<string>("data type");
408  if (!uinput->hasInputDataType(input_type))
409  {
410  std::cerr << "Input type:" + input_type + ", unavailable or misspelled."
411  << std::endl; // bad type
412  return adapter;
413  }
414 
415  vector<const zscalar_t *> weights;
416  std::vector<int> weightStrides;
417 
418  // get weights if any
419  if(uinput->hasUIWeights())
420  {
421  RCP<tMVector_t> vtx_weights = uinput->getUIWeights();
422  // copy to weight
423  size_t weightsPerRow = vtx_weights->getNumVectors();
424  for (size_t i = 0; i< weightsPerRow; i++) {
425  weights.push_back(vtx_weights->getData(i).getRawPtr());
426  weightStrides.push_back(1);
427  }
428  }
429 
430  // set adapter
431  if(input_type == "coordinates")
432  {
433  RCP<tMVector_t> data = uinput->getUICoordinates();
434  RCP<const tMVector_t> const_data = rcp_const_cast<const tMVector_t>(data);
435  if(weights.empty())
436  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(new Zoltan2::XpetraMultiVectorAdapter<tMVector_t>(const_data));
437  else
438  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(new Zoltan2::XpetraMultiVectorAdapter<tMVector_t>(const_data,weights,weightStrides));
439  }
440  else if(input_type == "tpetra_multivector")
441  {
442  int nvec = pList.get<int>("vector_dimension");
443  RCP<tMVector_t> data = uinput->getUITpetraMultiVector(nvec);
444  RCP<const tMVector_t> const_data = rcp_const_cast<const tMVector_t>(data);
445  if(weights.empty())
446  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(new Zoltan2::XpetraMultiVectorAdapter<tMVector_t>(const_data));
447  else
448  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(new Zoltan2::XpetraMultiVectorAdapter<tMVector_t>(const_data,weights,weightStrides));
449  }
450  else if(input_type == "xpetra_multivector")
451  {
452  int nvec = pList.get<int>("vector_dimension");
453  RCP<xMVector_t> data = uinput->getUIXpetraMultiVector(nvec);
454  RCP<const xMVector_t> const_data = rcp_const_cast<const xMVector_t>(data);
455  if(weights.empty())
456  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(new Zoltan2::XpetraMultiVectorAdapter<xMVector_t>(const_data));
457  else{
458  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(new Zoltan2::XpetraMultiVectorAdapter<xMVector_t>(const_data,weights,weightStrides));
459  }
460  }
461 #ifdef HAVE_EPETRA_DATA_TYPES
462 
463  else if(input_type == "epetra_multivector")
464  {
465  int nvec = pList.get<int>("vector_dimension");
466  RCP<Epetra_MultiVector> data = uinput->getUIEpetraMultiVector(nvec);
467  RCP<const Epetra_MultiVector> const_data = rcp_const_cast<const Epetra_MultiVector>(data);
468 
469  if(weights.empty())
470  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(
472  else
473  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(
474  new Zoltan2::XpetraMultiVectorAdapter<Epetra_MultiVector>(const_data,weights,weightStrides));
475  }
476 #endif
477 
478  if(adapter == nullptr)
479  std::cerr << "Input data chosen not compatible with xpetra multi-vector adapter." << std::endl;
480 
481  return adapter;
482 }
483 
484 
485 AdapterForTests::base_adapter_t * AdapterForTests::getXpetraCrsGraphAdapterForInput(
486  UserInputForTests *uinput,
487  const ParameterList &pList,
488  const RCP<const Comm<int> > &comm)
489 {
490 
491  AdapterForTests::base_adapter_t * adapter = nullptr;
492 
493  if(!pList.isParameter("data type"))
494  {
495  std::cerr << "Input data type unspecified" << std::endl;
496  return adapter;
497  }
498 
499  string input_type = pList.get<string>("data type");
500  if (!uinput->hasInputDataType(input_type))
501  {
502  std::cerr << "Input type: " + input_type + ", unavailable or misspelled."
503  << std::endl; // bad type
504  return adapter;
505  }
506 
507  vector<const zscalar_t *> vtx_weights;
508  vector<const zscalar_t *> edge_weights;
509  vector<int> vtx_weightStride;
510  vector<int> edge_weightStride;
511 
512  // get vtx weights if any
513  if(uinput->hasUIWeights())
514  {
515  RCP<tMVector_t> vtx_weights_tmp = uinput->getUIWeights();
516  // copy to weight
517  size_t weightsPerRow = vtx_weights_tmp->getNumVectors();
518  for (size_t i = 0; i< weightsPerRow; i++) {
519  vtx_weights.push_back(vtx_weights_tmp->getData(i).getRawPtr());
520  vtx_weightStride.push_back(1);
521  }
522  }
523 
524  // get edge weights if any
525  if(uinput->hasUIEdgeWeights())
526  {
527  RCP<tMVector_t> edge_weights_tmp = uinput->getUIEdgeWeights();
528  // copy to weight
529  size_t weightsPerRow = edge_weights_tmp->getNumVectors();
530  for (size_t i = 0; i< weightsPerRow; i++) {
531  edge_weights.push_back(edge_weights_tmp->getData(i).getRawPtr());
532  edge_weightStride.push_back(1);
533  }
534  }
535 
536 
537  // set adapter
538  if(input_type == "tpetra_crs_graph")
539  {
541 
542  RCP<tcrsGraph_t> data = uinput->getUITpetraCrsGraph();
543  RCP<const tcrsGraph_t> const_data = rcp_const_cast<const tcrsGraph_t>(data);
544  problem_t *ia = new problem_t(const_data,(int)vtx_weights.size(),(int)edge_weights.size());
545 
546  if(!vtx_weights.empty())
547  {
548  for(int i = 0; i < (int)vtx_weights.size(); i++)
549  ia->setVertexWeights(vtx_weights[i],vtx_weightStride[i],i);
550  }
551 
552  if(!edge_weights.empty())
553  {
554  for(int i = 0; i < (int)edge_weights.size(); i++)
555  ia->setEdgeWeights(edge_weights[i],edge_weightStride[i],i);
556  }
557 
558  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(ia);
559  }
560  else if(input_type == "xpetra_crs_graph")
561  {
563 
564  RCP<xcrsGraph_t> data = uinput->getUIXpetraCrsGraph();
565  RCP<const xcrsGraph_t> const_data = rcp_const_cast<const xcrsGraph_t>(data);
566  problem_t *ia = new problem_t(const_data, (int)vtx_weights.size(), (int)edge_weights.size());
567 
568  if(!vtx_weights.empty())
569  {
570  for(int i = 0; i < (int)vtx_weights.size(); i++)
571  ia->setVertexWeights(vtx_weights[i],vtx_weightStride[i],i);
572  }
573 
574  if(!edge_weights.empty())
575  {
576  for(int i = 0; i < (int)edge_weights.size(); i++)
577  ia->setEdgeWeights(edge_weights[i],edge_weightStride[i],i);
578  }
579 
580  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(ia);
581  }
582 #ifdef HAVE_EPETRA_DATA_TYPES
583 
584  else if(input_type == "epetra_crs_graph")
585  {
587 
588  RCP<Epetra_CrsGraph> data = uinput->getUIEpetraCrsGraph();
589  RCP<const Epetra_CrsGraph> const_data = rcp_const_cast<const Epetra_CrsGraph>(data);
590  problem_t *ia = new problem_t(const_data,(int)vtx_weights.size(),(int)edge_weights.size());
591 
592  if(!vtx_weights.empty())
593  {
594  for(int i = 0; i < (int)vtx_weights.size(); i++)
595  ia->setVertexWeights(vtx_weights[i],vtx_weightStride[i],i);
596  }
597 
598  if(!edge_weights.empty())
599  {
600  for(int i = 0; i < (int)edge_weights.size(); i++)
601  ia->setEdgeWeights(edge_weights[i],edge_weightStride[i],i);
602  }
603 
604  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(ia);
605 
606  }
607 #endif
608 
609  if(adapter == nullptr)
610  {
611  std::cerr << "Input data chosen not compatible with "
612  << "XpetraCrsGraph adapter." << std::endl;
613  return adapter;
614  }
615  else if (uinput->hasUICoordinates()) {
616  // make the coordinate adapter
617  // get an adapter for the coordinates
618  // need to make a copy of the plist and change the vector type
619  Teuchos::ParameterList pCopy(pList);
620  pCopy = pCopy.set<std::string>("data type","coordinates");
621 
622  AdapterForTests::base_adapter_t * ca = nullptr;
623  ca = getXpetraMVAdapterForInput(uinput,pCopy, comm);
624 
625  if(ca == nullptr)
626  {
627  std::cerr << "Failed to create coordinate vector adapter for "
628  << "XpetraCrsMatrix adapter." << std::endl;
629  return ca;
630  }
631 
632  // set the coordinate adapter
633  reinterpret_cast<AdapterForTests::xcrsGraph_adapter *>(adapter)->setCoordinateInput(reinterpret_cast<AdapterForTests::xpetra_mv_adapter *>(ca));
634  }
635  return adapter;
636 }
637 
638 
639 AdapterForTests::base_adapter_t * AdapterForTests::getXpetraCrsMatrixAdapterForInput(
640  UserInputForTests *uinput,
641  const ParameterList &pList,
642  const RCP<const Comm<int> > &comm)
643 {
644  AdapterForTests::base_adapter_t * adapter = nullptr;
645 
646  if(!pList.isParameter("data type"))
647  {
648  std::cerr << "Input data type unspecified" << std::endl;
649  return adapter;
650  }
651 
652  string input_type = pList.get<string>("data type");
653  if (!uinput->hasInputDataType(input_type))
654  {
655  std::cerr << "Input type:" + input_type + ", unavailable or misspelled."
656  << std::endl; // bad type
657  return adapter;
658  }
659 
660  vector<const zscalar_t *> weights;
661  vector<int> strides;
662 
663  // get weights if any
664  if(uinput->hasUIWeights())
665  {
666  if(comm->getRank() == 0) cout << "Have weights...." << endl;
667  RCP<tMVector_t> vtx_weights = uinput->getUIWeights();
668 
669  // copy to weight
670  int weightsPerRow = (int)vtx_weights->getNumVectors();
671  for (int i = 0; i< weightsPerRow; i++)
672  {
673  weights.push_back(vtx_weights->getData(i).getRawPtr());
674  strides.push_back(1);
675  }
676 
677  }
678 
679  // set adapter
680  if(input_type == "tpetra_crs_matrix")
681  {
682  if(comm->getRank() == 0) cout << "Make tpetra crs matrix adapter...." << endl;
683 
684  // get pointer to data
685  RCP<tcrsMatrix_t> data = uinput->getUITpetraCrsMatrix();
686  RCP<const tcrsMatrix_t> const_data = rcp_const_cast<const tcrsMatrix_t>(data); // const cast data
687 
688  // new adapter
689  xcrsMatrix_adapter *ia = new xcrsMatrix_adapter(const_data, (int)weights.size());
690 
691  // if we have weights set them
692  if(!weights.empty())
693  {
694  for(int i = 0; i < (int)weights.size(); i++)
695  ia->setWeights(weights[i],strides[i],i);
696  }
697 
698  // cast to base type
699  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(ia);
700 
701  }
702  else if(input_type == "xpetra_crs_matrix")
703  {
704  // type def this adapter type
706 
707  RCP<xcrsMatrix_t> data = uinput->getUIXpetraCrsMatrix();
708  RCP<const xcrsMatrix_t> const_data = rcp_const_cast<const xcrsMatrix_t>(data);
709 
710  // new adapter
711  problem_t *ia = new problem_t(const_data, (int)weights.size());
712 
713  // if we have weights set them
714  if(!weights.empty())
715  {
716  for(int i = 0; i < (int)weights.size(); i++)
717  ia->setWeights(weights[i],strides[i],i);
718  }
719 
720  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(ia);
721 
722  }
723 #ifdef HAVE_EPETRA_DATA_TYPES
724 
725  else if(input_type == "epetra_crs_matrix")
726  {
728 
729  RCP<Epetra_CrsMatrix> data = uinput->getUIEpetraCrsMatrix();
730  RCP<const Epetra_CrsMatrix> const_data = rcp_const_cast<const Epetra_CrsMatrix>(data);
731 
732  // new adapter
733  problem_t *ia = new problem_t(const_data, (int)weights.size());
734 
735  // if we have weights set them
736  if(!weights.empty())
737  {
738  for(int i = 0; i < (int)weights.size(); i++)
739  ia->setWeights(weights[i],strides[i],i);
740  }
741 
742  adapter = reinterpret_cast<AdapterForTests::base_adapter_t *>(ia);
743  }
744 #endif
745 
746  if(adapter == nullptr)
747  {
748  std::cerr << "Input data chosen not compatible with "
749  << "XpetraCrsMatrix adapter." << std::endl;
750  return adapter;
751  }
752  else if (uinput->hasUICoordinates()) {
753  // make the coordinate adapter
754  // get an adapter for the coordinates
755  // need to make a copy of the plist and change the vector type
756  Teuchos::ParameterList pCopy(pList);
757  pCopy = pCopy.set<std::string>("data type","coordinates");
758 
759  AdapterForTests::base_adapter_t * ca = nullptr;
760  ca = getXpetraMVAdapterForInput(uinput,pCopy,comm);
761 
762  if(ca == nullptr){
763  std::cerr << "Failed to create coordinate vector adapter for "
764  << "XpetraCrsMatrix adapter." << std::endl;
765  return ca;
766  }
767 
768  // set the coordinate adapter
769  reinterpret_cast<AdapterForTests::xcrsMatrix_adapter *>(adapter)->setCoordinateInput(reinterpret_cast<AdapterForTests::xpetra_mv_adapter *>(ca));
770  }
771  return adapter;
772 }
773 
774 
775 AdapterForTests::base_adapter_t * AdapterForTests::getBasicVectorAdapterForInput(
776  UserInputForTests *uinput,
777  const ParameterList &pList,
778  const RCP<const Comm<int> > &comm)
779 {
780 
781  AdapterForTests::basic_vector_adapter * ia = nullptr; // pointer for basic vector adapter
782 
783  if(!pList.isParameter("data type"))
784  {
785  std::cerr << "Input data type unspecified" << std::endl;
786  return nullptr;
787  }
788 
789  string input_type = pList.get<string>("data type");
790  if (!uinput->hasInputDataType(input_type))
791  {
792  std::cerr << "Input type:" + input_type + ", unavailable or misspelled."
793  << std::endl; // bad type
794  return nullptr;
795  }
796 
797  vector<const zscalar_t *> weights;
798  std::vector<int> weightStrides;
799  const zgno_t * globalIds;
800  zlno_t localCount = 0;
801 
802  // get weights if any
803  // get weights if any
804  if(uinput->hasUIWeights())
805  {
806  RCP<tMVector_t> vtx_weights = uinput->getUIWeights();
807  // copy to weight
808  size_t cols = vtx_weights->getNumVectors();
809  for (size_t i = 0; i< cols; i++) {
810  weights.push_back(vtx_weights->getData(i).getRawPtr());
811  weightStrides.push_back(1);
812  }
813  }
814 
815  // get vector stride
816  int stride = 1;
817  if(pList.isParameter("stride"))
818  stride = pList.get<int>("stride");
819 
820  if(input_type == "coordinates")
821  {
822  RCP<tMVector_t> data = uinput->getUICoordinates();
823  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
824  localCount = static_cast<zlno_t>(data->getLocalLength());
825 
826  // get strided data
827  vector<const zscalar_t *> coords;
828  vector<int> entry_strides;
829  AdapterForTests::InitializeVectorData(data,coords,entry_strides,stride);
830 
831  size_t dim = coords.size(); //BDD may need to add NULL for constructor call
832  size_t push_null = 3-dim;
833  for (size_t i = 0; i < push_null; i ++)
834  coords.push_back(NULL);
835 
836  if(weights.empty())
837  {
838  ia = new AdapterForTests::basic_vector_adapter(zlno_t(localCount),
839  globalIds,
840  coords[0],coords[1],coords[2],
841  stride, stride, stride);
842  }else{
843  ia = new AdapterForTests::basic_vector_adapter(zlno_t(localCount),
844  globalIds,
845  coords[0],coords[1],coords[2],
846  stride, stride, stride,
847  true,
848  weights[0],
849  weightStrides[0]);
850  }
851 
852  }
853  else if(input_type == "tpetra_vector")
854  {
855  RCP<tVector_t> data = uinput->getUITpetraVector();
856  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
857  localCount = static_cast<zlno_t>(data->getLocalLength());
858 
859  // get strided data
860  vector<const zscalar_t *> coords;
861  vector<int> entry_strides;
862  AdapterForTests::InitializeVectorData(data,coords,entry_strides,stride);
863 
864  if(weights.empty())
865  {
866  ia = new AdapterForTests::basic_vector_adapter(localCount, globalIds,
867  coords[0], entry_strides[0]);
868  }else{
869  ia = new AdapterForTests::basic_vector_adapter(localCount, globalIds,
870  coords[0], entry_strides[0],
871  true,
872  weights[0],
873  weightStrides[0]);
874 
875  }
876 
877  }
878  else if(input_type == "tpetra_multivector")
879  {
880  int nvec = pList.get<int>("vector_dimension");
881 
882  RCP<tMVector_t> data = uinput->getUITpetraMultiVector(nvec);
883  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
884  localCount = static_cast<zlno_t>(data->getLocalLength());
885 
886  // get strided data
887  vector<const zscalar_t *> coords;
888  vector<int> entry_strides;
889  AdapterForTests::InitializeVectorData(data,coords,entry_strides,stride);
890 
891  ia = new AdapterForTests::basic_vector_adapter(localCount, globalIds,
892  coords, entry_strides,
893  weights,weightStrides);
894 
895  }
896  else if(input_type == "xpetra_vector")
897  {
898  RCP<xVector_t> data = uinput->getUIXpetraVector();
899  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
900  localCount = static_cast<zlno_t>(data->getLocalLength());
901 
902  // get strided data
903  vector<const zscalar_t *> coords;
904  vector<int> entry_strides;
905  AdapterForTests::InitializeVectorData(data,coords,entry_strides,stride);
906 
907  if(weights.empty())
908  {
909  ia = new AdapterForTests::basic_vector_adapter(localCount, globalIds,
910  coords[0], entry_strides[0]);
911  }else{
912  ia = new AdapterForTests::basic_vector_adapter(localCount, globalIds,
913  coords[0], entry_strides[0],
914  true,
915  weights[0],
916  weightStrides[0]);
917 
918  }
919  }
920  else if(input_type == "xpetra_multivector")
921  {
922  int nvec = pList.get<int>("vector_dimension");
923  RCP<xMVector_t> data = uinput->getUIXpetraMultiVector(nvec);
924  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
925  localCount = static_cast<zlno_t>(data->getLocalLength());
926 
927  // get strided data
928  vector<const zscalar_t *> coords;
929  vector<int> entry_strides;
930  AdapterForTests::InitializeVectorData(data,coords,entry_strides,stride);
931  if(comm->getRank() == 0) cout << "size of entry strides: " << entry_strides.size() << endl;
932  if(comm->getRank() == 0) cout << "size of coords: " << coords.size() << endl;
933 
934  // make vector!
935  ia = new AdapterForTests::basic_vector_adapter(localCount, globalIds,
936  coords, entry_strides,
937  weights,weightStrides);
938  }
939 
940 #ifdef HAVE_EPETRA_DATA_TYPES
941  else if(input_type == "epetra_vector")
942  {
943  RCP<Epetra_Vector> data = uinput->getUIEpetraVector();
944  globalIds = (zgno_t *)data->Map().MyGlobalElements();
945  localCount = static_cast<zlno_t>(data->MyLength());
946 
947  // get strided data
948  vector<const zscalar_t *> coords;
949  vector<int> entry_strides;
950  AdapterForTests::InitializeEpetraVectorData(data,coords,entry_strides,stride);
951 
952  if(weights.empty())
953  {
954  ia = new AdapterForTests::basic_vector_adapter(localCount, globalIds,
955  coords[0], entry_strides[0]);
956  }else{
957  ia = new AdapterForTests::basic_vector_adapter(localCount, globalIds,
958  coords[0], entry_strides[0],
959  true,
960  weights[0],
961  weightStrides[0]);
962 
963  }
964 
965  // delete [] epetravectors;
966  }
967  else if(input_type == "epetra_multivector")
968  {
969  int nvec = pList.get<int>("vector_dimension");
970  RCP<Epetra_MultiVector> data = uinput->getUIEpetraMultiVector(nvec);
971  globalIds = (zgno_t *)data->Map().MyGlobalElements();
972  localCount = data->MyLength();
973 
974  vector<const zscalar_t *> coords;
975  vector<int> entry_strides;
976  AdapterForTests::InitializeEpetraVectorData(data,coords,entry_strides,stride);
977 
978  // make vector!
979  ia = new AdapterForTests::basic_vector_adapter(localCount, globalIds,
980  coords, entry_strides,
981  weights,weightStrides);
982 
983  }
984 
985 #endif
986 
987  if(localCount == 0){
988  if(ia != nullptr) delete ia;
989  return nullptr;
990  }
991  return reinterpret_cast<AdapterForTests::base_adapter_t *>(ia);
992 
993 }
994 
995 template <typename T>
996 void AdapterForTests::InitializeVectorData(const RCP<T> &data,
997  vector<const zscalar_t *> &coords,
998  vector<int> & strides,
999  int stride)
1000 {
1001  // set up adapter data
1002  size_t localCount = data->getLocalLength();
1003  size_t nvecs = data->getNumVectors();
1004  size_t vecsize = data->getNumVectors() * data->getLocalLength();
1005 // printf("Number of vectors by data: %zu\n", nvecs);
1006  // printf("Size of data: %zu\n", vecsize);
1007 
1008  ArrayRCP<zscalar_t> *petravectors =
1009  new ArrayRCP<zscalar_t>[nvecs];
1010 
1011  // printf("Getting t-petra vectors...\n");
1012  for (size_t i = 0; i < nvecs; i++)
1013  petravectors[i] = data->getDataNonConst(i);
1014 
1015  // debugging
1016  // for (size_t i = 0; i < nvecs; i++){
1017  // printf("Tpetra vector %zu: {",i);
1018  //
1019  // for (size_t j = 0; j < localCount; j++)
1020  // {
1021  // printf("%1.2g ",petravectors[i][j]);
1022  // }
1023  // printf("}\n");
1024  // }
1025 
1026  size_t idx = 0;
1027  zscalar_t *coordarr = new zscalar_t[vecsize];
1028 
1029  if(stride == 1 || stride != (int)nvecs)
1030  {
1031  for (size_t i = 0; i < nvecs; i++) {
1032  for (size_t j = 0; j < localCount; j++) {
1033  coordarr[idx++] = petravectors[i][j];
1034  }
1035  }
1036  }else
1037  {
1038  for (size_t j = 0; j < localCount; j++) {
1039  for (size_t i = 0; i < nvecs; i++) {
1040  coordarr[idx++] = petravectors[i][j];
1041  }
1042  }
1043  }
1044 
1045  // debugging
1046  // printf("Made coordarr : {");
1047  // for (zlno_t i = 0; i < vecsize; i++){
1048  // printf("%1.2g ",coordarr[i]);
1049  // }
1050  // printf("}\n");
1051 
1052  // always build for dim 3
1053  coords = std::vector<const zscalar_t *>(nvecs);
1054  strides = std::vector<int>(nvecs);
1055 
1056  for (size_t i = 0; i < nvecs; i++) {
1057  if(stride == 1)
1058  coords[i] = &coordarr[i*localCount];
1059  else
1060  coords[i] = &coordarr[i];
1061 
1062  strides[i] = stride;
1063  }
1064 
1065  // debugging
1066  // printf("Made coords...\n");
1067  // for (size_t i = 0; i < nvecs; i++){
1068  // const zscalar_t * tmp = coords[i];
1069  // printf("coord %zu: {",i);
1070  // for(size_t j = 0; j < localCount; j++)
1071  // {
1072  // printf("%1.2g ", tmp[j]);
1073  // }
1074  // printf("}\n");
1075  // }
1076 
1077  // printf("clean up coordarr and tpetravectors...\n\n\n");
1078  delete [] petravectors;
1079 }
1080 
1081 #ifdef HAVE_EPETRA_DATA_TYPES
1082 
1083 template <typename T>
1084 void AdapterForTests::InitializeEpetraVectorData(const RCP<T> &data,
1085  vector<const zscalar_t *> &coords,
1086  vector<int> & strides,
1087  int stride){
1088  size_t localCount = data->MyLength();
1089  size_t nvecs = data->NumVectors();
1090  size_t vecsize = nvecs * localCount;
1091 
1092  // printf("Number of vectors by data: %zu\n", nvecs);
1093  // printf("Size of data: %zu\n", vecsize);
1094 
1095  vector<zscalar_t *> epetravectors(nvecs);
1096  zscalar_t ** arr;
1097  // printf("get data from epetra vector..\n");
1098  data->ExtractView(&arr);
1099 
1100  for(size_t k = 0; k < nvecs; k++)
1101  {
1102  epetravectors[k] = arr[k];
1103  }
1104 
1105  size_t idx = 0;
1107 
1108  if(stride == 1 || stride != (int)nvecs)
1109  {
1110  for (size_t i = 0; i < nvecs; i++) {
1111  for (size_t j = 0; j < localCount; j++) {
1112  coordarr[idx++] = epetravectors[i][j];
1113  }
1114  }
1115  }else
1116  {
1117  for (size_t j = 0; j < localCount; j++) {
1118  for (size_t i = 0; i < nvecs; i++) {
1119  coordarr[idx++] = epetravectors[i][j];
1120  }
1121  }
1122  }
1123 
1124  // debugging
1125 // printf("Made coordarr : {");
1126 // for (zlno_t i = 0; i < vecsize; i++){
1127 // printf("%1.2g ",coordarr[i]);
1128 // }
1129 // printf("}\n");
1130 
1131  coords = std::vector<const zscalar_t *>(nvecs);
1132  strides = std::vector<int>(nvecs);
1133 
1134  for (size_t i = 0; i < nvecs; i++) {
1135  if(stride == 1)
1136  coords[i] = &coordarr[i*localCount];
1137  else
1138  coords[i] = &coordarr[i];
1139 
1140  strides[i] = stride;
1141  }
1142 
1143 // printf("Made coords...\n");
1144 // for (size_t i = 0; i < nvecs; i++){
1145 // const zscalar_t * tmp = coords[i];
1146 // printf("coord %zu: {",i);
1147 // for(size_t j = 0; j < localCount; j++)
1148 // {
1149 // printf("%1.2g ", tmp[j]);
1150 // }
1151 // printf("}\n");
1152 // }
1153 
1154 }
1155 #endif
1156 
1157 
1158 // pamgen adapter
1160 AdapterForTests::getPamgenMeshAdapterForInput(UserInputForTests *uinput,
1161  const ParameterList &pList,
1162  const RCP<const Comm<int> > &comm)
1163 {
1164 #ifdef HAVE_ZOLTAN2_PAMGEN
1165  pamgen_adapter_t * ia = nullptr; // pointer for basic vector adapter
1166  if(uinput->hasPamgenMesh())
1167  {
1168 
1169  if(uinput->hasPamgenMesh())
1170  {
1171 // if(comm->getRank() == 0) cout << "Have pamgen mesh, constructing adapter...." << endl;
1172  ia = new pamgen_adapter_t(*(comm.get()), "region");
1173 // if(comm->getRank() == 0)
1174 // ia->print(0);
1175  }
1176  }else{
1177  std::cerr << "Pamgen mesh is unavailable for PamgenMeshAdapter!"
1178  << std::endl;
1179  }
1180 
1181  return reinterpret_cast<AdapterForTests::base_adapter_t *>(ia);
1182 #else
1183  throw std::runtime_error("Pamgen input requested but Trilinos is not "
1184  "built with Pamgen");
1185 #endif
1186 }
1187 #endif
1188 
1189 
InputTraits< User >::scalar_t scalar_t
Generate input for testing purposes.
UserInputForTests::tMVector_t tMVector_t
Zoltan2::XpetraMultiVectorAdapter< tMVector_t > xpetra_mv_adapter
Zoltan2::XpetraCrsGraphAdapter< tcrsGraph_t, tMVector_t > xcrsGraph_adapter
Zoltan2::BaseAdapter< userTypes_t > base_adapter_t
Xpetra::CrsGraph< zlno_t, zgno_t, znode_t > xcrsGraph_t
Tpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > tcrsMatrix_t
Defines Parameter related enumerators, declares functions.
Provides access for Zoltan2 to Xpetra::CrsMatrix data.
UserInputForTests::xMVector_t xMVector_t
double zscalar_t
A simple class that can be the User template argument for an InputAdapter.
UserInputForTests::tcrsGraph_t tcrsGraph_t
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tMVector_t
RCP< tMVector_t > getUIWeights()
Provides access for Zoltan2 to Xpetra::CrsGraph data.
static ArrayRCP< ArrayRCP< zscalar_t > > weights
Defines the PamgenMeshAdapter class.
int zlno_t
UserInputForTests::tcrsMatrix_t tcrsMatrix_t
Zoltan2::XpetraCrsMatrixAdapter< tcrsMatrix_t, tMVector_t > xcrsMatrix_adapter
UserInputForTests::xcrsGraph_t xcrsGraph_t
Zoltan2::BasicIdentifierAdapter< userTypes_t > basic_id_t
Defines the XpetraMultiVectorAdapter.
Zoltan2::BasicVectorAdapter< tMVector_t > basic_vector_adapter
RCP< xcrsGraph_t > getUIXpetraCrsGraph()
Defines XpetraCrsGraphAdapter class.
RCP< xVector_t > getUIXpetraVector()
This class represents a collection of global Identifiers and their associated weights, if any.
RCP< xMVector_t > getUIXpetraMultiVector(int nvec)
Defines the XpetraCrsMatrixAdapter class.
bool hasInputDataType(const string &input_type)
RCP< tMVector_t > getUITpetraMultiVector(int nvec)
RCP< tVector_t > getUITpetraVector()
Defines the EvaluatePartition class.
RCP< xcrsMatrix_t > getUIXpetraCrsMatrix()
Xpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > xcrsMatrix_t
RCP< tcrsGraph_t > getUITpetraCrsGraph()
Zoltan2::BasicVectorAdapter< tMVector_t > pamgen_adapter_t
Xpetra::Vector< zscalar_t, zlno_t, zgno_t, znode_t > xVector_t
UserInputForTests::xVector_t xVector_t
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
Xpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > xMVector_t
An adapter for Xpetra::MultiVector.
RCP< tcrsMatrix_t > getUITpetraCrsMatrix()
RCP< tMVector_t > getUIEdgeWeights()
UserInputForTests::xcrsMatrix_t xcrsMatrix_t
int zgno_t
Defines the BasicIdentifierAdapter class.
BaseAdapter defines methods required by all Adapters.
Tpetra::Vector< zscalar_t, zlno_t, zgno_t, znode_t > tVector_t
void setWeights(const scalar_t *weightVal, int stride, int idx=0)
Specify a weight for each entity of the primaryEntityType.
static base_adapter_t * getAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP< const Comm< int > > &comm)
A class method for constructing an input adapter defind in a parameter list.
Defines the PartitioningProblem class.
RCP< tMVector_t > getUICoordinates()
Zoltan2::BasicUserTypes< zscalar_t, zlno_t, zgno_t > userTypes_t
Defines the BasicVectorAdapter class.
This class represents a mesh.
UserInputForTests::tVector_t tVector_t
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > tcrsGraph_t