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 using namespace Zoltan2_TestingFramework;
93 
94 // helper struct to store both an adapter and the coordinate adapter
96 {
98  Zoltan2::VectorAdapter<tMVector_t> * coordinateAdapter = nullptr; // may be NULL if coordinates are not available
99 };
100 
101 
102 /* \brief A class for constructing Zoltan2 input adapters */
104 public:
105 
106 #ifdef HAVE_ZOLTAN2_PAMGEN
108 #else
109  // This typedef exists only to satisfy the compiler.
110  // PamgenMeshAdapter cannot be used when Trilinos is not built with Pamgen
112 #endif
113 
123  static AdapterWithOptionalCoordinateAdapter getAdapterForInput(
124  UserInputForTests *uinput, const ParameterList &pList,
125  const RCP<const Comm<int> > &comm);
126 
127 private:
136  static base_adapter_t*
137  getBasicIdentiferAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
138 
147  static base_adapter_t*
148  getXpetraMVAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
149 
159  getXpetraCrsGraphAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
160 
170  getXpetraCrsMatrixAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
171 
180  static base_adapter_t*
181  getBasicVectorAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
182 
191  static base_adapter_t*
192  getPamgenMeshAdapterForInput(UserInputForTests *uinput, const ParameterList &pList, const RCP<const Comm<int> > &comm);
193 
194 
203  template <typename T>
204  static void InitializeVectorData(const RCP<T> &data,
205  vector<const zscalar_t *> &coords,
206  vector<int> & strides,
207  int stride);
208 
209 #ifdef HAVE_EPETRA_DATA_TYPES
210 
218  template <typename T>
219  static void InitializeEpetraVectorData(const RCP<T> &data,
220  vector<const zscalar_t *> &coords,
221  vector<int> & strides,
222  int stride);
223 #endif
224 };
225 
226 
228  UserInputForTests *uinput,
229  const ParameterList &pList,
230  const RCP<const Comm<int> > &comm)
231 {
232  AdapterWithOptionalCoordinateAdapter adapters; // input adapter
233 
234  if(!pList.isParameter("input adapter"))
235  {
236  std::cerr << "Input adapter unspecified" << std::endl;
237  return adapters;
238  }
239 
240 
241  // pick method for chosen adapter
242  string adapter_name = pList.get<string>("input adapter");
243 
244  if(adapter_name == "BasicIdentifier")
245  adapters.mainAdapter = AdapterForTests::getBasicIdentiferAdapterForInput(uinput, pList, comm);
246  else if(adapter_name == "XpetraMultiVector")
247  adapters.mainAdapter = AdapterForTests::getXpetraMVAdapterForInput(uinput, pList, comm);
248  else if(adapter_name == "XpetraCrsGraph")
249  adapters = getXpetraCrsGraphAdapterForInput(uinput,pList, comm);
250  else if(adapter_name == "XpetraCrsMatrix")
251  adapters = getXpetraCrsMatrixAdapterForInput(uinput,pList, comm);
252  else if(adapter_name == "BasicVector")
253  adapters.mainAdapter = getBasicVectorAdapterForInput(uinput,pList, comm);
254  else if(adapter_name == "PamgenMesh")
255  adapters.mainAdapter = getPamgenMeshAdapterForInput(uinput,pList, comm);
256  else
257  std::cerr << "Input adapter type: " << adapter_name
258  << ", is unavailable, or misspelled." << std::endl;
259 
260  return adapters;
261 }
262 
263 
264 Zoltan2_TestingFramework::base_adapter_t * AdapterForTests::getBasicIdentiferAdapterForInput(UserInputForTests *uinput,
265  const ParameterList &pList,
266  const RCP<const Comm<int> > &comm)
267 {
268 
269  if(!pList.isParameter("data type"))
270  {
271  std::cerr << "Input data type unspecified" << std::endl;
272  return nullptr;
273  }
274 
275  string input_type = pList.get<string>("data type"); // get the input type
276 
277  if (!uinput->hasInputDataType(input_type))
278  {
279  std::cerr << "Input type: " + input_type + " unavailable or misspelled."
280  << std::endl; // bad type
281  return nullptr;
282  }
283 
284  vector<const zscalar_t *> weights;
285  std::vector<int> weightStrides;
286  const zgno_t *globalIds = NULL;
287  size_t localCount = 0;
288 
289  // get weights if any
290  // get weights if any
291  if(uinput->hasUIWeights())
292  {
293  RCP<tMVector_t> vtx_weights = uinput->getUIWeights();
294  // copy to weight
295  size_t cols = vtx_weights->getNumVectors();
296  for (size_t i = 0; i< cols; i++) {
297  weights.push_back(vtx_weights->getData(i).getRawPtr());
298  weightStrides.push_back((int)vtx_weights->getStride());
299  }
300  }
301 
302  if(input_type == "coordinates")
303  {
304  RCP<tMVector_t> data = uinput->getUICoordinates();
305  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
306  localCount = data->getLocalLength();
307  }
308  else if(input_type == "tpetra_vector")
309  {
310  RCP<tVector_t> data = uinput->getUITpetraVector();
311  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
312  localCount = data->getLocalLength();
313  }
314  else if(input_type == "tpetra_multivector")
315  {
316  int nvec = pList.get<int>("vector_dimension");
317  RCP<tMVector_t> data = uinput->getUITpetraMultiVector(nvec);
318  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
319  localCount = data->getLocalLength();
320  }
321  else if(input_type == "tpetra_crs_graph")
322  {
323  RCP<tcrsGraph_t> data = uinput->getUITpetraCrsGraph();
324  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
325  localCount = data->getNodeNumCols();
326  }
327  else if(input_type == "tpetra_crs_matrix")
328  {
329  RCP<tcrsMatrix_t> data = uinput->getUITpetraCrsMatrix();
330  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
331  localCount = data->getNodeNumCols();
332  }
333  else if(input_type == "xpetra_vector")
334  {
335  RCP<xVector_t> data = uinput->getUIXpetraVector();
336  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
337  localCount = data->getLocalLength();
338  }
339  else if(input_type == "xpetra_multivector")
340  {
341  int nvec = pList.get<int>("vector_dimension");
342  RCP<xMVector_t> data = uinput->getUIXpetraMultiVector(nvec);
343  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
344  localCount = data->getLocalLength();
345  }
346  else if(input_type == "xpetra_crs_graph")
347  {
348  RCP<xcrsGraph_t> data = uinput->getUIXpetraCrsGraph();
349  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
350  localCount = data->getNodeNumCols();
351  }
352  else if(input_type == "xpetra_crs_matrix")
353  {
354  RCP<xcrsMatrix_t> data = uinput->getUIXpetraCrsMatrix();
355  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
356  localCount = data->getNodeNumCols();
357  }
358 #ifdef HAVE_EPETRA_DATA_TYPES
359  else if(input_type == "epetra_vector")
360  {
361  RCP<Epetra_Vector> data = uinput->getUIEpetraVector();
362  globalIds = (zgno_t *)data->Map().MyGlobalElements();
363  localCount = data->MyLength();
364  }
365  else if(input_type == "epetra_multivector")
366  {
367  int nvec = pList.get<int>("vector_dimension");
368  RCP<Epetra_MultiVector> data = uinput->getUIEpetraMultiVector(nvec);
369  globalIds = (zgno_t *)data->Map().MyGlobalElements();
370  localCount = data->MyLength();
371  }
372  else if(input_type == "epetra_crs_graph")
373  {
374  RCP<Epetra_CrsGraph> data = uinput->getUIEpetraCrsGraph();
375  globalIds = (zgno_t *)data->Map().MyGlobalElements();
376  localCount = data->NumMyCols();
377  }
378  else if(input_type == "epetra_crs_matrix")
379  {
380  RCP<Epetra_CrsMatrix> data = uinput->getUIEpetraCrsMatrix();
381  globalIds = (zgno_t *)data->Map().MyGlobalElements();
382  localCount = data->NumMyCols();
383  }
384 #endif
385 
386  return reinterpret_cast<Zoltan2_TestingFramework::base_adapter_t *>(
388  zlno_t(localCount),
389  globalIds,
390  weights,weightStrides));
391 }
392 
393 
394 Zoltan2_TestingFramework::base_adapter_t * AdapterForTests::getXpetraMVAdapterForInput(
395  UserInputForTests *uinput,
396  const ParameterList &pList,
397  const RCP<const Comm<int> > &comm)
398 {
399  Zoltan2_TestingFramework::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<Zoltan2_TestingFramework::base_adapter_t *>(new Zoltan2::XpetraMultiVectorAdapter<tMVector_t>(const_data));
437  else
438  adapter = reinterpret_cast<Zoltan2_TestingFramework::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<Zoltan2_TestingFramework::base_adapter_t *>(new Zoltan2::XpetraMultiVectorAdapter<tMVector_t>(const_data));
447  else
448  adapter = reinterpret_cast<Zoltan2_TestingFramework::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<Zoltan2_TestingFramework::base_adapter_t *>(new Zoltan2::XpetraMultiVectorAdapter<xMVector_t>(const_data));
457  else{
458  adapter = reinterpret_cast<Zoltan2_TestingFramework::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<Zoltan2_TestingFramework::base_adapter_t *>(
472  else
473  adapter = reinterpret_cast<Zoltan2_TestingFramework::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 AdapterWithOptionalCoordinateAdapter AdapterForTests::getXpetraCrsGraphAdapterForInput(
486  UserInputForTests *uinput,
487  const ParameterList &pList,
488  const RCP<const Comm<int> > &comm)
489 {
490 
492 
493  if(!pList.isParameter("data type"))
494  {
495  std::cerr << "Input data type unspecified" << std::endl;
496  return adapters;
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 adapters;
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  adapter_t *ia = new adapter_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  adapters.mainAdapter = reinterpret_cast<Zoltan2_TestingFramework::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  adapter_t *ia = new adapter_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  adapters.mainAdapter = reinterpret_cast<Zoltan2_TestingFramework::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  adapter_t *ia = new adapter_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  adapters.mainAdapter = reinterpret_cast<Zoltan2_TestingFramework::base_adapter_t *>(ia);
605  }
606 #endif
607 
608  if(adapters.mainAdapter == nullptr)
609  {
610  std::cerr << "Input data chosen not compatible with "
611  << "XpetraCrsGraph adapter." << std::endl;
612  return adapters;
613  }
614  else if (uinput->hasUICoordinates()) {
615  // make the coordinate adapter
616  // get an adapter for the coordinates
617  // need to make a copy of the plist and change the vector type
618  Teuchos::ParameterList pCopy(pList);
619  pCopy = pCopy.set<std::string>("data type","coordinates");
620 
622  ca = getXpetraMVAdapterForInput(uinput,pCopy, comm);
623 
624  if(ca == nullptr)
625  {
626  std::cerr << "Failed to create coordinate vector adapter for "
627  << "XpetraCrsMatrix adapter." << std::endl;
628  return adapters;
629  }
630 
631  adapters.coordinateAdapter = reinterpret_cast<Zoltan2_TestingFramework::xpetra_mv_adapter *>(ca);
632 
633  // set the coordinate adapter
635  (adapters.mainAdapter)->setCoordinateInput(
636  adapters.coordinateAdapter);
637 
638  }
639  return adapters;
640 }
641 
642 
643 AdapterWithOptionalCoordinateAdapter AdapterForTests::getXpetraCrsMatrixAdapterForInput(
644  UserInputForTests *uinput,
645  const ParameterList &pList,
646  const RCP<const Comm<int> > &comm)
647 {
649 
650  if(!pList.isParameter("data type"))
651  {
652  std::cerr << "Input data type unspecified" << std::endl;
653  return adapters;
654  }
655 
656  string input_type = pList.get<string>("data type");
657  if (!uinput->hasInputDataType(input_type))
658  {
659  std::cerr << "Input type:" + input_type + ", unavailable or misspelled."
660  << std::endl; // bad type
661  return adapters;
662  }
663 
664  vector<const zscalar_t *> weights;
665  vector<int> strides;
666 
667  // get weights if any
668  if(uinput->hasUIWeights())
669  {
670  if(comm->getRank() == 0) cout << "Have weights...." << endl;
671  RCP<tMVector_t> vtx_weights = uinput->getUIWeights();
672 
673  // copy to weight
674  int weightsPerRow = (int)vtx_weights->getNumVectors();
675  for (int i = 0; i< weightsPerRow; i++)
676  {
677  weights.push_back(vtx_weights->getData(i).getRawPtr());
678  strides.push_back(1);
679  }
680 
681  }
682 
683  // set adapter
684  if(input_type == "tpetra_crs_matrix")
685  {
686  if(comm->getRank() == 0) cout << "Make tpetra crs matrix adapter...." << endl;
687 
688  // get pointer to data
689  RCP<tcrsMatrix_t> data = uinput->getUITpetraCrsMatrix();
690  RCP<const tcrsMatrix_t> const_data = rcp_const_cast<const tcrsMatrix_t>(data); // const cast data
691 
692  // new adapter
693  xcrsMatrix_adapter *ia = new xcrsMatrix_adapter(const_data, (int)weights.size());
694 
695  // if we have weights set them
696  if(!weights.empty())
697  {
698  for(int i = 0; i < (int)weights.size(); i++)
699  ia->setWeights(weights[i],strides[i],i);
700  }
701 
702  // cast to base type
703  adapters.mainAdapter = reinterpret_cast<Zoltan2_TestingFramework::base_adapter_t *>(ia);
704 
705  }
706  else if(input_type == "xpetra_crs_matrix")
707  {
708  // type def this adapter type
710 
711  RCP<xcrsMatrix_t> data = uinput->getUIXpetraCrsMatrix();
712  RCP<const xcrsMatrix_t> const_data = rcp_const_cast<const xcrsMatrix_t>(data);
713 
714  // new adapter
715  adapter_t *ia = new adapter_t(const_data, (int)weights.size());
716 
717  // if we have weights set them
718  if(!weights.empty())
719  {
720  for(int i = 0; i < (int)weights.size(); i++)
721  ia->setWeights(weights[i],strides[i],i);
722  }
723 
724  adapters.mainAdapter = reinterpret_cast<Zoltan2_TestingFramework::base_adapter_t *>(ia);
725 
726  }
727 #ifdef HAVE_EPETRA_DATA_TYPES
728 
729  else if(input_type == "epetra_crs_matrix")
730  {
732 
733  RCP<Epetra_CrsMatrix> data = uinput->getUIEpetraCrsMatrix();
734  RCP<const Epetra_CrsMatrix> const_data = rcp_const_cast<const Epetra_CrsMatrix>(data);
735 
736  // new adapter
737  adapter_t *ia = new adapter_t(const_data, (int)weights.size());
738 
739  // if we have weights set them
740  if(!weights.empty())
741  {
742  for(int i = 0; i < (int)weights.size(); i++)
743  ia->setWeights(weights[i],strides[i],i);
744  }
745 
746  adapters.mainAdapter = reinterpret_cast<Zoltan2_TestingFramework::base_adapter_t *>(ia);
747  }
748 #endif
749 
750  if(adapters.mainAdapter == nullptr)
751  {
752  std::cerr << "Input data chosen not compatible with "
753  << "XpetraCrsMatrix adapter." << std::endl;
754  return adapters;
755  }
756  else if (uinput->hasUICoordinates()) {
757  // make the coordinate adapter
758  // get an adapter for the coordinates
759  // need to make a copy of the plist and change the vector type
760  Teuchos::ParameterList pCopy(pList);
761  pCopy = pCopy.set<std::string>("data type","coordinates");
762 
764  ca = getXpetraMVAdapterForInput(uinput,pCopy,comm);
765 
766  if(ca == nullptr){
767  std::cerr << "Failed to create coordinate vector adapter for "
768  << "XpetraCrsMatrix adapter." << std::endl;
769  return adapters;
770  }
771 
772  // set the coordinate adapter
773  adapters.coordinateAdapter =
774  reinterpret_cast<Zoltan2_TestingFramework::xpetra_mv_adapter *>(ca);
775 
777  (adapters.mainAdapter)->setCoordinateInput(adapters.coordinateAdapter);
778 
779 
780  }
781  return adapters;
782 }
783 
784 
785 Zoltan2_TestingFramework::base_adapter_t * AdapterForTests::getBasicVectorAdapterForInput(
786  UserInputForTests *uinput,
787  const ParameterList &pList,
788  const RCP<const Comm<int> > &comm)
789 {
790 
792 
793  if(!pList.isParameter("data type"))
794  {
795  std::cerr << "Input data type unspecified" << std::endl;
796  return nullptr;
797  }
798 
799  string input_type = pList.get<string>("data type");
800  if (!uinput->hasInputDataType(input_type))
801  {
802  std::cerr << "Input type:" + input_type + ", unavailable or misspelled."
803  << std::endl; // bad type
804  return nullptr;
805  }
806 
807  std::vector<const zscalar_t *> weights;
808  std::vector<int> weightStrides;
809  const zgno_t * globalIds;
810  zlno_t localCount = 0;
811 
812  // get weights if any
813  // get weights if any
814  if(uinput->hasUIWeights())
815  {
816  RCP<tMVector_t> vtx_weights = uinput->getUIWeights();
817  // copy to weight
818  size_t cols = vtx_weights->getNumVectors();
819  for (size_t i = 0; i< cols; i++) {
820  weights.push_back(vtx_weights->getData(i).getRawPtr());
821  weightStrides.push_back(1);
822  }
823  }
824 
825  // get vector stride
826  int stride = 1;
827  if(pList.isParameter("stride"))
828  stride = pList.get<int>("stride");
829 
830  if(input_type == "coordinates")
831  {
832  RCP<tMVector_t> data = uinput->getUICoordinates();
833  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
834  localCount = static_cast<zlno_t>(data->getLocalLength());
835 
836  // get strided data
837  std::vector<const zscalar_t *> coords;
838  std::vector<int> entry_strides;
839  AdapterForTests::InitializeVectorData(data,coords,entry_strides,stride);
840 
841 
842  if (weights.empty()) {
843  size_t dim = coords.size(); //BDD add NULL for constructor call
844  size_t push_null = 3-dim;
845  for (size_t i = 0; i < push_null; i ++) coords.push_back(NULL);
847  zlno_t(localCount),
848  globalIds,
849  coords[0],
850  coords[1],coords[2],
851  stride, stride, stride);
852  } else if (weights.size() == 1) {
853  size_t dim = coords.size(); //BDD add NULL for constructor call
854  size_t push_null = 3-dim;
855  for (size_t i = 0; i < push_null; i ++) coords.push_back(NULL);
857  zlno_t(localCount),
858  globalIds,
859  coords[0],
860  coords[1],coords[2],
861  stride, stride, stride,
862  true,
863  weights[0],
864  weightStrides[0]);
865  } else { // More than one weight per ID
867  zlno_t(localCount),
868  globalIds,
869  coords, entry_strides,
870  weights, weightStrides);
871  }
872  }
873  else if(input_type == "tpetra_vector")
874  {
875  RCP<tVector_t> data = uinput->getUITpetraVector();
876  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
877  localCount = static_cast<zlno_t>(data->getLocalLength());
878 
879  // get strided data
880  vector<const zscalar_t *> coords;
881  vector<int> entry_strides;
882  AdapterForTests::InitializeVectorData(data,coords,entry_strides,stride);
883 
884  if(weights.empty())
885  {
886  ia = new Zoltan2_TestingFramework::basic_vector_adapter(localCount, globalIds,
887  coords[0], entry_strides[0]);
888  }else{
889  ia = new Zoltan2_TestingFramework::basic_vector_adapter(localCount, globalIds,
890  coords[0], entry_strides[0],
891  true,
892  weights[0],
893  weightStrides[0]);
894 
895  }
896 
897  }
898  else if(input_type == "tpetra_multivector")
899  {
900  int nvec = pList.get<int>("vector_dimension");
901 
902  RCP<tMVector_t> data = uinput->getUITpetraMultiVector(nvec);
903  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
904  localCount = static_cast<zlno_t>(data->getLocalLength());
905 
906  // get strided data
907  vector<const zscalar_t *> coords;
908  vector<int> entry_strides;
909  AdapterForTests::InitializeVectorData(data,coords,entry_strides,stride);
910 
911  ia = new Zoltan2_TestingFramework::basic_vector_adapter(localCount, globalIds,
912  coords, entry_strides,
913  weights,weightStrides);
914 
915  }
916  else if(input_type == "xpetra_vector")
917  {
918  RCP<xVector_t> data = uinput->getUIXpetraVector();
919  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
920  localCount = static_cast<zlno_t>(data->getLocalLength());
921 
922  // get strided data
923  vector<const zscalar_t *> coords;
924  vector<int> entry_strides;
925  AdapterForTests::InitializeVectorData(data,coords,entry_strides,stride);
926 
927  if(weights.empty())
928  {
929  ia = new Zoltan2_TestingFramework::basic_vector_adapter(localCount, globalIds,
930  coords[0], entry_strides[0]);
931  }else{
932  ia = new Zoltan2_TestingFramework::basic_vector_adapter(localCount, globalIds,
933  coords[0], entry_strides[0],
934  true,
935  weights[0],
936  weightStrides[0]);
937 
938  }
939  }
940  else if(input_type == "xpetra_multivector")
941  {
942  int nvec = pList.get<int>("vector_dimension");
943  RCP<xMVector_t> data = uinput->getUIXpetraMultiVector(nvec);
944  globalIds = (zgno_t *)data->getMap()->getNodeElementList().getRawPtr();
945  localCount = static_cast<zlno_t>(data->getLocalLength());
946 
947  // get strided data
948  vector<const zscalar_t *> coords;
949  vector<int> entry_strides;
950  AdapterForTests::InitializeVectorData(data,coords,entry_strides,stride);
951  if(comm->getRank() == 0) cout << "size of entry strides: " << entry_strides.size() << endl;
952  if(comm->getRank() == 0) cout << "size of coords: " << coords.size() << endl;
953 
954  // make vector!
955  ia = new Zoltan2_TestingFramework::basic_vector_adapter(localCount, globalIds,
956  coords, entry_strides,
957  weights,weightStrides);
958  }
959 
960 #ifdef HAVE_EPETRA_DATA_TYPES
961  else if(input_type == "epetra_vector")
962  {
963  RCP<Epetra_Vector> data = uinput->getUIEpetraVector();
964  globalIds = (zgno_t *)data->Map().MyGlobalElements();
965  localCount = static_cast<zlno_t>(data->MyLength());
966 
967  // get strided data
968  vector<const zscalar_t *> coords;
969  vector<int> entry_strides;
970  AdapterForTests::InitializeEpetraVectorData(data,coords,entry_strides,stride);
971 
972  if(weights.empty())
973  {
974  ia = new Zoltan2_TestingFramework::basic_vector_adapter(localCount, globalIds,
975  coords[0], entry_strides[0]);
976  }else{
977  ia = new Zoltan2_TestingFramework::basic_vector_adapter(localCount, globalIds,
978  coords[0], entry_strides[0],
979  true,
980  weights[0],
981  weightStrides[0]);
982 
983  }
984 
985  // delete [] epetravectors;
986  }
987  else if(input_type == "epetra_multivector")
988  {
989  int nvec = pList.get<int>("vector_dimension");
990  RCP<Epetra_MultiVector> data = uinput->getUIEpetraMultiVector(nvec);
991  globalIds = (zgno_t *)data->Map().MyGlobalElements();
992  localCount = data->MyLength();
993 
994  vector<const zscalar_t *> coords;
995  vector<int> entry_strides;
996  AdapterForTests::InitializeEpetraVectorData(data,coords,entry_strides,stride);
997 
998  // make vector!
999  ia = new Zoltan2_TestingFramework::basic_vector_adapter(localCount, globalIds,
1000  coords, entry_strides,
1001  weights,weightStrides);
1002 
1003  }
1004 
1005 #endif
1006 
1007  return reinterpret_cast<Zoltan2_TestingFramework::base_adapter_t *>(ia);
1008 
1009 }
1010 
1011 template <typename T>
1012 void AdapterForTests::InitializeVectorData(const RCP<T> &data,
1013  vector<const zscalar_t *> &coords,
1014  vector<int> & strides,
1015  int stride)
1016 {
1017  // set up adapter data
1018  size_t localCount = data->getLocalLength();
1019  size_t nvecs = data->getNumVectors();
1020  size_t vecsize = data->getNumVectors() * data->getLocalLength();
1021 // printf("Number of vectors by data: %zu\n", nvecs);
1022  // printf("Size of data: %zu\n", vecsize);
1023 
1024  ArrayRCP<zscalar_t> *petravectors =
1025  new ArrayRCP<zscalar_t>[nvecs];
1026 
1027  // printf("Getting t-petra vectors...\n");
1028  for (size_t i = 0; i < nvecs; i++)
1029  petravectors[i] = data->getDataNonConst(i);
1030 
1031  // debugging
1032  // for (size_t i = 0; i < nvecs; i++){
1033  // printf("Tpetra vector %zu: {",i);
1034  //
1035  // for (size_t j = 0; j < localCount; j++)
1036  // {
1037  // printf("%1.2g ",petravectors[i][j]);
1038  // }
1039  // printf("}\n");
1040  // }
1041 
1042  size_t idx = 0;
1043  zscalar_t *coordarr = new zscalar_t[vecsize];
1044 
1045  if(stride == 1 || stride != (int)nvecs)
1046  {
1047  for (size_t i = 0; i < nvecs; i++) {
1048  for (size_t j = 0; j < localCount; j++) {
1049  coordarr[idx++] = petravectors[i][j];
1050  }
1051  }
1052  }else
1053  {
1054  for (size_t j = 0; j < localCount; j++) {
1055  for (size_t i = 0; i < nvecs; i++) {
1056  coordarr[idx++] = petravectors[i][j];
1057  }
1058  }
1059  }
1060 
1061  // debugging
1062  // printf("Made coordarr : {");
1063  // for (zlno_t i = 0; i < vecsize; i++){
1064  // printf("%1.2g ",coordarr[i]);
1065  // }
1066  // printf("}\n");
1067 
1068  // always build for dim 3
1069  coords = std::vector<const zscalar_t *>(nvecs);
1070  strides = std::vector<int>(nvecs);
1071 
1072  for (size_t i = 0; i < nvecs; i++) {
1073  if(stride == 1)
1074  coords[i] = &coordarr[i*localCount];
1075  else
1076  coords[i] = &coordarr[i];
1077 
1078  strides[i] = stride;
1079  }
1080 
1081  // debugging
1082  // printf("Made coords...\n");
1083  // for (size_t i = 0; i < nvecs; i++){
1084  // const zscalar_t * tmp = coords[i];
1085  // printf("coord %zu: {",i);
1086  // for(size_t j = 0; j < localCount; j++)
1087  // {
1088  // printf("%1.2g ", tmp[j]);
1089  // }
1090  // printf("}\n");
1091  // }
1092 
1093  // printf("clean up coordarr and tpetravectors...\n\n\n");
1094  delete [] petravectors;
1095 }
1096 
1097 #ifdef HAVE_EPETRA_DATA_TYPES
1098 
1099 template <typename T>
1100 void AdapterForTests::InitializeEpetraVectorData(const RCP<T> &data,
1101  vector<const zscalar_t *> &coords,
1102  vector<int> & strides,
1103  int stride){
1104  size_t localCount = data->MyLength();
1105  size_t nvecs = data->NumVectors();
1106  size_t vecsize = nvecs * localCount;
1107 
1108  // printf("Number of vectors by data: %zu\n", nvecs);
1109  // printf("Size of data: %zu\n", vecsize);
1110 
1111  vector<zscalar_t *> epetravectors(nvecs);
1112  zscalar_t ** arr;
1113  // printf("get data from epetra vector..\n");
1114  data->ExtractView(&arr);
1115 
1116  for(size_t k = 0; k < nvecs; k++)
1117  {
1118  epetravectors[k] = arr[k];
1119  }
1120 
1121  size_t idx = 0;
1123 
1124  if(stride == 1 || stride != (int)nvecs)
1125  {
1126  for (size_t i = 0; i < nvecs; i++) {
1127  for (size_t j = 0; j < localCount; j++) {
1128  coordarr[idx++] = epetravectors[i][j];
1129  }
1130  }
1131  }else
1132  {
1133  for (size_t j = 0; j < localCount; j++) {
1134  for (size_t i = 0; i < nvecs; i++) {
1135  coordarr[idx++] = epetravectors[i][j];
1136  }
1137  }
1138  }
1139 
1140  // debugging
1141 // printf("Made coordarr : {");
1142 // for (zlno_t i = 0; i < vecsize; i++){
1143 // printf("%1.2g ",coordarr[i]);
1144 // }
1145 // printf("}\n");
1146 
1147  coords = std::vector<const zscalar_t *>(nvecs);
1148  strides = std::vector<int>(nvecs);
1149 
1150  for (size_t i = 0; i < nvecs; i++) {
1151  if(stride == 1)
1152  coords[i] = &coordarr[i*localCount];
1153  else
1154  coords[i] = &coordarr[i];
1155 
1156  strides[i] = stride;
1157  }
1158 
1159 // printf("Made coords...\n");
1160 // for (size_t i = 0; i < nvecs; i++){
1161 // const zscalar_t * tmp = coords[i];
1162 // printf("coord %zu: {",i);
1163 // for(size_t j = 0; j < localCount; j++)
1164 // {
1165 // printf("%1.2g ", tmp[j]);
1166 // }
1167 // printf("}\n");
1168 // }
1169 
1170 }
1171 #endif
1172 
1173 
1174 // pamgen adapter
1176 AdapterForTests::getPamgenMeshAdapterForInput(UserInputForTests *uinput,
1177  const ParameterList &pList,
1178  const RCP<const Comm<int> > &comm)
1179 {
1180 #ifdef HAVE_ZOLTAN2_PAMGEN
1181  pamgen_adapter_t * ia = nullptr; // pointer for basic vector adapter
1182  if(uinput->hasPamgenMesh())
1183  {
1184 
1185  if(uinput->hasPamgenMesh())
1186  {
1187 // if(comm->getRank() == 0) cout << "Have pamgen mesh, constructing adapter...." << endl;
1188  ia = new pamgen_adapter_t(*(comm.get()), "region");
1189 // if(comm->getRank() == 0)
1190 // ia->print(0);
1191  }
1192  }else{
1193  std::cerr << "Pamgen mesh is unavailable for PamgenMeshAdapter!"
1194  << std::endl;
1195  }
1196 
1197  return reinterpret_cast<Zoltan2_TestingFramework::base_adapter_t *>(ia);
1198 #else
1199  throw std::runtime_error("Pamgen input requested but Trilinos is not "
1200  "built with Pamgen");
1201 #endif
1202 }
1203 #endif
1204 
1205 
InputTraits< User >::scalar_t scalar_t
RCP< tMVector_t > getUITpetraMultiVector(int nvec)
Tpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > tcrsMatrix_t
Definition: GraphModel.cpp:84
Generate input for testing purposes.
RCP< tcrsMatrix_t > getUITpetraCrsMatrix()
Defines Parameter related enumerators, declares functions.
Provides access for Zoltan2 to Xpetra::CrsMatrix data.
RCP< tVector_t > getUITpetraVector()
double zscalar_t
static AdapterWithOptionalCoordinateAdapter 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.
Provides access for Zoltan2 to Xpetra::CrsGraph data.
static ArrayRCP< ArrayRCP< zscalar_t > > weights
Defines the PamgenMeshAdapter class.
int zlno_t
RCP< tMVector_t > getUIWeights()
Zoltan2::BasicVectorAdapter< tMVector_t > basic_vector_adapter
Defines the XpetraMultiVectorAdapter.
RCP< xcrsMatrix_t > getUIXpetraCrsMatrix()
Defines XpetraCrsGraphAdapter class.
RCP< xVector_t > getUIXpetraVector()
Defines the XpetraCrsMatrixAdapter class.
bool hasInputDataType(const string &input_type)
Xpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > xMVector_t
RCP< tMVector_t > getUIEdgeWeights()
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > tcrsGraph_t
Definition: GraphModel.cpp:85
Defines the EvaluatePartition class.
Xpetra::CrsGraph< zlno_t, zgno_t, znode_t > xcrsGraph_t
Zoltan2::BasicVectorAdapter< tMVector_t > pamgen_adapter_t
Xpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > xcrsMatrix_t
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
An adapter for Xpetra::MultiVector.
int zgno_t
Defines the BasicIdentifierAdapter class.
RCP< xMVector_t > getUIXpetraMultiVector(int nvec)
BaseAdapter defines methods required by all Adapters.
Zoltan2_TestingFramework::base_adapter_t * mainAdapter
Tpetra::MultiVector< double, int, int > tMVector_t
RCP< tMVector_t > getUICoordinates()
Zoltan2::XpetraCrsMatrixAdapter< tcrsMatrix_t, tMVector_t > xcrsMatrix_adapter
Defines the PartitioningProblem class.
RCP< xcrsGraph_t > getUIXpetraCrsGraph()
Defines the BasicVectorAdapter class.
This class represents a mesh.
Zoltan2::VectorAdapter< tMVector_t > * coordinateAdapter
Zoltan2::BasicIdentifierAdapter< userTypes_t > basic_id_t
Zoltan2::BasicVectorAdapter< tMVector_t > pamgen_adapter_t
RCP< tcrsGraph_t > getUITpetraCrsGraph()