50 #ifndef USERINPUTFORTESTS
51 #define USERINPUTFORTESTS
57 #include <Tpetra_MultiVector.hpp>
58 #include <Tpetra_CrsMatrix.hpp>
59 #include <Tpetra_Map.hpp>
60 #include <Xpetra_Vector.hpp>
61 #include <Xpetra_CrsMatrix.hpp>
62 #include <Xpetra_CrsGraph.hpp>
64 #include <MatrixMarket_Tpetra.hpp>
66 #ifdef HAVE_ZOLTAN2_GALERI
67 #include <Galeri_XpetraProblemFactory.hpp>
68 #include <Galeri_XpetraParameters.hpp>
71 #include <Tpetra_KokkosCompat_DefaultNode.hpp>
77 #include <TpetraExt_MatrixMatrix_def.hpp>
84 using Teuchos::ArrayRCP;
85 using Teuchos::ArrayView;
90 using Teuchos::rcp_const_cast;
91 using Teuchos::ParameterList;
92 using namespace Zoltan2_TestingFramework;
130 typedef Tpetra::Map<zlno_t, zgno_t, znode_t>
map_t;
131 typedef Tpetra::Export<zlno_t, zgno_t, znode_t>
export_t;
132 typedef Tpetra::Import<zlno_t, zgno_t, znode_t>
import_t;
156 const RCP<
const Comm<int> > &c,
bool debugInfo=
false,
157 bool distributeInput=
true);
176 const RCP<
const Comm<int> > &c,
bool debugInfo=
false,
177 bool distributeInput=
true);
195 const RCP<
const Comm<int> > &c);
199 static void getUIRandomData(
unsigned int seed,
zlno_t length,
202 RCP<tMVector_t> getUICoordinates();
204 RCP<tMVector_t> getUIWeights();
206 RCP<tMVector_t> getUIEdgeWeights();
208 RCP<tcrsMatrix_t> getUITpetraCrsMatrix();
210 RCP<tcrsGraph_t> getUITpetraCrsGraph();
212 RCP<tVector_t> getUITpetraVector();
214 RCP<tMVector_t> getUITpetraMultiVector(
int nvec);
216 RCP<xcrsMatrix_t> getUIXpetraCrsMatrix();
218 RCP<xcrsGraph_t> getUIXpetraCrsGraph();
220 RCP<xVector_t> getUIXpetraVector();
222 RCP<xMVector_t> getUIXpetraMultiVector(
int nvec);
224 #ifdef HAVE_ZOLTAN2_PAMGEN
225 PamgenMesh * getPamGenMesh(){
return this->pamgen_mesh.operator->();}
228 #ifdef HAVE_EPETRA_DATA_TYPES
229 RCP<Epetra_CrsGraph> getUIEpetraCrsGraph();
231 RCP<Epetra_CrsMatrix> getUIEpetraCrsMatrix();
233 RCP<Epetra_Vector> getUIEpetraVector();
235 RCP<Epetra_MultiVector> getUIEpetraMultiVector(
int nvec);
239 bool hasInputDataType(
const string &input_type);
241 bool hasUICoordinates();
245 bool hasUIEdgeWeights();
247 bool hasUITpetraCrsMatrix();
249 bool hasUITpetraCrsGraph();
251 bool hasUITpetraVector();
253 bool hasUITpetraMultiVector();
255 bool hasUIXpetraCrsMatrix();
257 bool hasUIXpetraCrsGraph();
259 bool hasUIXpetraVector();
261 bool hasUIXpetraMultiVector();
263 bool hasPamgenMesh();
264 #ifdef HAVE_EPETRA_DATA_TYPES
265 bool hasUIEpetraCrsGraph();
267 bool hasUIEpetraCrsMatrix();
269 bool hasUIEpetraVector();
271 bool hasUIEpetraMultiVector();
279 const RCP<const Comm<int> > tcomm_;
282 #ifdef HAVE_ZOLTAN2_PAMGEN
283 RCP<PamgenMesh> pamgen_mesh;
286 RCP<tcrsMatrix_t> M_;
287 RCP<xcrsMatrix_t> xM_;
289 RCP<tMVector_t> xyz_;
290 RCP<tMVector_t> vtxWeights_;
291 RCP<tMVector_t> edgWeights_;
293 #ifdef HAVE_EPETRA_DATA_TYPES
294 RCP<const Epetra_Comm> ecomm_;
295 RCP<Epetra_CrsMatrix> eM_;
296 RCP<Epetra_CrsGraph> eG_;
304 void readMatrixMarketFile(
string path,
string testData,
bool distributeInput =
true);
310 bool distributeInput);
316 void readZoltanTestData(
string path,
string testData,
317 bool distributeInput);
320 RCP<tcrsMatrix_t> modifyMatrixGIDs(RCP<tcrsMatrix_t> &in);
321 inline zgno_t newID(
const zgno_t id) {
return id * 2 + 10001; }
324 void getUIChacoGraph(FILE *fptr,
bool haveAssign, FILE *assignFile,
325 string name,
bool distributeInput);
328 void getUIChacoCoords(FILE *fptr,
string name);
334 static const int CHACO_LINE_LENGTH=200;
335 char chaco_line[CHACO_LINE_LENGTH];
340 double chaco_read_val(FILE* infile,
int *end_flag);
341 int chaco_read_int(FILE* infile,
int *end_flag);
342 void chaco_flush_line(FILE*);
345 int chaco_input_graph(FILE *fin,
const char *inname,
int **start,
346 int **adjacency,
int *nvtxs,
int *nVwgts,
347 float **vweights,
int *nEwgts,
float **eweights);
350 int chaco_input_geom(FILE *fingeom,
const char *geomname,
int nvtxs,
351 int *igeom,
double **x,
double **y,
double **z);
354 int chaco_input_assign(FILE *finassign,
const char *assignname,
int nvtxs,
362 void readGeometricGenTestData(
string path,
string testData);
366 ParameterList &geoparams);
371 const string& delimiters =
" \f\n\r\t\v" );
374 const string& delimiters =
" \f\n\r\t\v" );
377 const string& delimiters =
" \f\n\r\t\v" );
381 void readPamgenMeshFile(
string path,
string testData);
382 #ifdef HAVE_ZOLTAN2_PAMGEN
383 void setPamgenAdjacencyGraph();
384 void setPamgenCoordinateMV();
389 const RCP<
const Comm<int> > &c,
390 bool debugInfo,
bool distributeInput):
391 verbose_(debugInfo), tcomm_(c), havePamgenMesh(false),
392 M_(), xM_(), xyz_(), vtxWeights_(), edgWeights_(),
393 #ifdef HAVE_EPETRA_DATA_TYPES
394 ecomm_(), eM_(), eG_(),
396 chaco_offset(0), chaco_break_pnt(CHACO_LINE_LENGTH)
398 bool zoltan1 =
false;
399 string::size_type loc = path.find(
"/zoltan/test/");
400 if (loc != string::npos)
404 readZoltanTestData(path, testData, distributeInput);
406 readMatrixMarketFile(path, testData);
408 #ifdef HAVE_EPETRA_DATA_TYPES
409 ecomm_ = Xpetra::toEpetra(c);
415 const RCP<
const Comm<int> > &c,
417 bool distributeInput):
418 verbose_(debugInfo), tcomm_(c), havePamgenMesh(false),
419 M_(), xM_(), xyz_(), vtxWeights_(), edgWeights_(),
420 #ifdef HAVE_EPETRA_DATA_TYPES
421 ecomm_(), eM_(), eG_(),
423 chaco_offset(0), chaco_break_pnt(CHACO_LINE_LENGTH)
425 if (matrixType.size() == 0){
431 matrixType = string(
"Laplace1D");
433 matrixType = string(
"Laplace2D");
435 matrixType = string(
"Laplace3D");
437 throw std::runtime_error(
"input");
439 if (verbose_ && tcomm_->getRank() == 0)
440 std::cout <<
"UserInputForTests, Matrix type : " << matrixType << std::endl;
443 buildCrsMatrix(x, y, z, matrixType, distributeInput);
445 #ifdef HAVE_EPETRA_DATA_TYPES
446 ecomm_ = Xpetra::toEpetra(c);
451 const RCP<
const Comm<int> > &c):
452 tcomm_(c), havePamgenMesh(false),
453 M_(), xM_(), xyz_(), vtxWeights_(), edgWeights_(),
454 #ifdef HAVE_EPETRA_DATA_TYPES
455 ecomm_(), eM_(), eG_(),
457 chaco_offset(0), chaco_break_pnt(CHACO_LINE_LENGTH)
461 bool distributeInput =
true, debugInfo =
true;
463 if(pList.isParameter(
"distribute input"))
464 distributeInput = pList.get<
bool>(
"distribute input");
466 if(pList.isParameter(
"debug"))
467 debugInfo = pList.get<
bool>(
"debug");
468 this->verbose_ = debugInfo;
470 if(pList.isParameter(
"input file"))
475 if(pList.isParameter(
"input path"))
476 path = pList.get<
string>(
"input path");
478 string testData = pList.get<
string>(
"input file");
484 if(pList.isParameter(
"file type") && pList.get<
string>(
"file type") ==
"Geometric Generator")
486 else if(pList.isParameter(
"file type") && pList.get<
string>(
"file type") ==
"Pamgen")
490 else if(pList.isParameter(
"file type") && pList.get<
string>(
"file type") ==
"Chaco")
494 switch (file_format) {
495 case GEOMGEN: readGeometricGenTestData(path,testData);
break;
496 case PAMGEN: readPamgenMeshFile(path,testData);
break;
497 case CHACO: readZoltanTestData(path, testData, distributeInput);
break;
498 default: readMatrixMarketFile(path, testData, distributeInput);
break;
501 }
else if(pList.isParameter(
"x") || pList.isParameter(
"y") || pList.isParameter(
"z")){
505 if(pList.isParameter(
"x")) x = pList.get<
int>(
"x");
506 if(pList.isParameter(
"y")) y = pList.get<
int>(
"y");
507 if(pList.isParameter(
"z")) z = pList.get<
int>(
"z");
509 string problemType =
"";
510 if(pList.isParameter(
"equation type")) problemType = pList.get<
string>(
"equation type");
512 if (problemType.size() == 0){
518 problemType = string(
"Laplace1D");
520 problemType = string(
"Laplace2D");
522 problemType = string(
"Laplace3D");
524 throw std::runtime_error(
"input");
526 if (verbose_ && tcomm_->getRank() == 0)
527 std::cout <<
"UserInputForTests, Matrix type : " << problemType << std::endl;
531 buildCrsMatrix(x, y, z, problemType, distributeInput);
534 std::cerr <<
"Input file block undefined!" << std::endl;
537 #ifdef HAVE_EPETRA_DATA_TYPES
538 ecomm_ = Xpetra::toEpetra(c);
547 throw std::runtime_error(
"could not read coord file");
564 throw std::runtime_error(
"could not read mtx file");
571 throw std::runtime_error(
"could not read mtx file");
572 return rcp_const_cast<
tcrsGraph_t>(M_->getCrsGraph());
577 RCP<tVector_t> V = rcp(
new tVector_t(M_->getRowMap(), 1));
585 RCP<tMVector_t> mV = rcp(
new tMVector_t(M_->getRowMap(), nvec));
594 throw std::runtime_error(
"could not read mtx file");
601 throw std::runtime_error(
"could not read mtx file");
602 return rcp_const_cast<
xcrsGraph_t>(xM_->getCrsGraph());
616 #ifdef HAVE_EPETRA_DATA_TYPES
617 RCP<Epetra_CrsGraph> UserInputForTests::getUIEpetraCrsGraph()
620 throw std::runtime_error(
"could not read mtx file");
621 RCP<const tcrsGraph_t> tgraph = M_->getCrsGraph();
622 RCP<const Tpetra::Map<zlno_t, zgno_t> > trowMap = tgraph->getRowMap();
623 RCP<const Tpetra::Map<zlno_t, zgno_t> > tcolMap = tgraph->getColMap();
625 int nElts =
static_cast<int>(trowMap->getGlobalNumElements());
626 int nMyElts =
static_cast<int>(trowMap->getLocalNumElements());
628 ArrayView<const int> gids = trowMap->getLocalElementList();
630 Epetra_BlockMap erowMap(nElts, nMyElts,
631 gids.getRawPtr(), 1, base, *ecomm_);
633 Array<int> rowSize(nMyElts);
634 for (
int i=0; i < nMyElts; i++){
635 rowSize[i] =
static_cast<int>(M_->getNumEntriesInLocalRow(i));
638 size_t maxRow = M_->getLocalMaxNumRowEntries();
639 Array<int> colGids(maxRow);
641 typename tcrsGraph_t::local_inds_host_view_type colLid;
642 eG_ = rcp(
new Epetra_CrsGraph(Copy, erowMap,
643 rowSize.getRawPtr(),
true));
645 for (
int i=0; i < nMyElts; i++){
646 tgraph->getLocalRowView(i, colLid);
647 for (
size_t j=0; j < colLid.extent(0); j++)
648 colGids[j] = tcolMap->getGlobalElement(colLid[j]);
649 eG_->InsertGlobalIndices(gids[i], rowSize[i], colGids.getRawPtr());
655 RCP<Epetra_CrsMatrix> UserInputForTests::getUIEpetraCrsMatrix()
658 throw std::runtime_error(
"could not read mtx file");
659 RCP<Epetra_CrsGraph> egraph = getUIEpetraCrsGraph();
660 eM_ = rcp(
new Epetra_CrsMatrix(Copy, *egraph));
662 size_t maxRow = M_->getLocalMaxNumRowEntries();
663 int nrows = egraph->NumMyRows();
664 const Epetra_BlockMap &rowMap = egraph->RowMap();
665 const Epetra_BlockMap &colMap = egraph->ColMap();
666 Array<int> colGid(maxRow);
667 for (
int i=0; i < nrows; i++){
668 typename tcrsMatrix_t::local_inds_host_view_type colLid;
669 typename tcrsMatrix_t::values_host_view_type nz;
670 M_->getLocalRowView(i, colLid, nz);
671 size_t rowSize = colLid.size();
672 int rowGid = rowMap.GID(i);
673 for (
size_t j=0; j < rowSize; j++){
674 colGid[j] = colMap.GID(colLid[j]);
676 eM_->InsertGlobalValues(rowGid, (
int)rowSize, nz.data(), colGid.getRawPtr());
682 RCP<Epetra_Vector> UserInputForTests::getUIEpetraVector()
684 RCP<Epetra_CrsGraph> egraph = getUIEpetraCrsGraph();
685 RCP<Epetra_Vector> V = rcp(
new Epetra_Vector(egraph->RowMap()));
690 RCP<Epetra_MultiVector> UserInputForTests::getUIEpetraMultiVector(
int nvec)
692 RCP<Epetra_CrsGraph> egraph = getUIEpetraCrsGraph();
693 RCP<Epetra_MultiVector> mV =
694 rcp(
new Epetra_MultiVector(egraph->RowMap(), nvec));
704 this->hasUITpetraCrsMatrix() || \
705 this->hasUITpetraCrsGraph() || \
706 this->hasPamgenMesh();
711 if(input_type ==
"coordinates")
713 else if(input_type ==
"tpetra_vector")
715 else if(input_type ==
"tpetra_multivector")
717 else if(input_type ==
"tpetra_crs_graph")
719 else if(input_type ==
"tpetra_crs_matrix")
721 else if(input_type ==
"xpetra_vector")
723 else if(input_type ==
"xpetra_multivector")
725 else if(input_type ==
"xpetra_crs_graph")
727 else if(input_type ==
"xpetra_crs_matrix")
729 #ifdef HAVE_EPETRA_DATA_TYPES
730 else if(input_type ==
"epetra_vector")
731 return this->hasUIEpetraVector();
732 else if(input_type ==
"epetra_multivector")
733 return this->hasUIEpetraMultiVector();
734 else if(input_type ==
"epetra_crs_graph")
735 return this->hasUIEpetraCrsGraph();
736 else if(input_type ==
"epetra_crs_matrix")
737 return this->hasUIEpetraCrsMatrix();
745 return xyz_.is_null() ?
false :
true;
750 return vtxWeights_.is_null() ?
false :
true;
755 return edgWeights_.is_null() ?
false :
true;
760 return M_.is_null() ?
false :
true;
765 return M_.is_null() ?
false :
true;
780 return M_.is_null() ?
false :
true;
785 return M_.is_null() ?
false :
true;
800 return this->havePamgenMesh;
803 #ifdef HAVE_EPETRA_DATA_TYPES
804 bool UserInputForTests::hasUIEpetraCrsGraph()
806 return M_.is_null() ?
false :
true;
809 bool UserInputForTests::hasUIEpetraCrsMatrix()
811 return hasUIEpetraCrsGraph();
814 bool UserInputForTests::hasUIEpetraVector()
816 return hasUIEpetraCrsGraph();
819 bool UserInputForTests::hasUIEpetraMultiVector()
821 return hasUIEpetraCrsGraph();
827 ArrayView<ArrayRCP<zscalar_t > > data)
832 size_t dim = data.size();
833 for (
size_t i=0; i < dim; i++){
836 throw (std::bad_alloc());
837 data[i] = Teuchos::arcp(tmp, 0, length,
true);
840 zscalar_t scalingFactor = (max-min) / RAND_MAX;
842 for (
size_t i=0; i < dim; i++){
844 for (
zlno_t j=0; j < length; j++)
845 *x++ = min + (
zscalar_t(rand()) * scalingFactor);
851 string UserInputForTests::trim_right_copy(
853 const string& delimiters)
855 return s.substr( 0, s.find_last_not_of( delimiters ) + 1 );
858 string UserInputForTests::trim_left_copy(
860 const string& delimiters)
862 return s.substr( s.find_first_not_of( delimiters ) );
865 string UserInputForTests::trim_copy(
867 const string& delimiters)
869 return trim_left_copy( trim_right_copy( s, delimiters ), delimiters );
872 void UserInputForTests::readGeometricGenTestData(
string path,
876 std::ostringstream
fname;
877 fname << path <<
"/" << testData <<
".txt";
879 if (verbose_ && tcomm_->getRank() == 0)
880 std::cout <<
"UserInputForTests, Read: " << fname.str() << std::endl;
882 Teuchos::ParameterList geoparams(
"geo params");
883 readGeoGenParams(fname.str(),geoparams);
888 int coord_dim = gg->getCoordinateDimension();
889 int numWeightsPerCoord = gg->getNumWeights();
890 zlno_t numLocalPoints = gg->getNumLocalCoords();
891 zgno_t numGlobalPoints = gg->getNumGlobalCoords();
895 for(
int i = 0; i < coord_dim; ++i){
896 coords[i] =
new zscalar_t[numLocalPoints];
900 gg->getLocalCoordinatesCopy(coords);
904 if (numWeightsPerCoord) {
906 weight =
new zscalar_t * [numWeightsPerCoord];
907 for(
int i = 0; i < numWeightsPerCoord; ++i){
908 weight[i] =
new zscalar_t[numLocalPoints];
912 gg->getLocalWeightsCopy(weight);
919 RCP<Tpetra::Map<zlno_t, zgno_t, znode_t> > mp =
920 rcp(
new Tpetra::Map<zlno_t, zgno_t, znode_t>(numGlobalPoints, numLocalPoints, 0, this->tcomm_));
923 Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > coordView(coord_dim);
924 for (
int i=0; i < coord_dim; i++){
925 if(numLocalPoints > 0){
926 Teuchos::ArrayView<const zscalar_t> a(coords[i], numLocalPoints);
930 Teuchos::ArrayView<const zscalar_t> a;
936 xyz_ = RCP<tMVector_t>(
new
941 if (numWeightsPerCoord) {
943 Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > weightView(numWeightsPerCoord);
944 for (
int i=0; i < numWeightsPerCoord; i++){
945 if(numLocalPoints > 0){
946 Teuchos::ArrayView<const zscalar_t> a(weight[i], numLocalPoints);
950 Teuchos::ArrayView<const zscalar_t> a;
955 vtxWeights_ = RCP<tMVector_t>(
new tMVector_t(mp, weightView.view(0, numWeightsPerCoord),
956 numWeightsPerCoord));
960 void UserInputForTests::readGeoGenParams(
string paramFileName,
961 ParameterList &geoparams){
965 std::string input =
"";
967 for(
int i = 0; i < 25000; ++i){
972 if(this->tcomm_->getRank() == 0){
974 std::fstream inParam(paramFileName.c_str());
981 std::string tmp =
"";
982 getline (inParam,tmp);
983 while (!inParam.eof()){
985 tmp = trim_copy(tmp);
990 getline (inParam,tmp);
993 for (
size_t i = 0; i < input.size(); ++i){
1001 int size = (int)input.size();
1005 this->tcomm_->broadcast(0,
sizeof(
int), (
char*) &size);
1007 throw "File " + paramFileName +
" cannot be opened.";
1009 this->tcomm_->broadcast(0, size, inp);
1010 std::istringstream inParam(inp);
1012 getline (inParam,str);
1013 while (!inParam.eof()){
1014 if(str[0] != param_comment){
1015 size_t pos = str.find(
'=');
1016 if(pos == string::npos){
1017 throw "Invalid Line:" + str +
" in parameter file";
1019 string paramname = trim_copy(str.substr(0,pos));
1020 string paramvalue = trim_copy(str.substr(pos + 1));
1021 geoparams.set(paramname, paramvalue);
1023 getline (inParam,str);
1028 RCP<tcrsMatrix_t> UserInputForTests::modifyMatrixGIDs(
1029 RCP<tcrsMatrix_t> &inMatrix
1038 RCP<const map_t> inMap = inMatrix->getRowMap();
1040 size_t nRows = inMap->getLocalNumElements();
1041 auto inRows = inMap->getMyGlobalIndices();
1042 Teuchos::Array<zgno_t> outRows(nRows);
1043 for (
size_t i = 0; i < nRows; i++) {
1044 outRows[i] = newID(inRows[i]);
1048 RCP<map_t> outMap = rcp(
new map_t(nGlobalRows, outRows(), 0,
1051 #ifdef INCLUDE_LENGTHY_OUTPUT
1054 std::cout << inMap->getComm()->getRank() <<
" KDDKDD "
1055 <<
"nGlobal " << inMap->getGlobalNumElements() <<
" "
1056 << outMap->getGlobalNumElements() <<
"; "
1057 <<
"nLocal " << inMap->getLocalNumElements() <<
" "
1058 << outMap->getLocalNumElements() <<
"; "
1060 std::cout << inMap->getComm()->getRank() <<
" KDDKDD ";
1061 for (
size_t i = 0; i < nRows; i++)
1062 std::cout <<
"(" << inMap->getMyGlobalIndices()[i] <<
", "
1063 << outMap->getMyGlobalIndices()[i] <<
") ";
1064 std::cout << std::endl;
1066 #endif // INCLUDE_LENGTHY_OUTPUT
1070 size_t rowLen = inMatrix->getLocalMaxNumRowEntries();
1071 RCP<tcrsMatrix_t> outMatrix = rcp(
new tcrsMatrix_t(outMap, rowLen));
1073 typename tcrsMatrix_t::nonconst_global_inds_host_view_type indices(
"Indices", rowLen);
1074 typename tcrsMatrix_t::nonconst_values_host_view_type values(
"Values", rowLen);
1076 for (
size_t i = 0; i < nRows; i++) {
1078 zgno_t inGid = inMap->getGlobalElement(i);
1079 inMatrix->getGlobalRowCopy(inGid, indices, values, nEntries);
1080 for (
size_t j = 0; j < nEntries; j++)
1081 indices[j] = newID(indices[j]);
1083 zgno_t outGid = outMap->getGlobalElement(i);
1084 ArrayView<zgno_t> Indices(indices.data(), nEntries);
1085 ArrayView<zscalar_t> Values(values.data(), nEntries);
1086 outMatrix->insertGlobalValues(outGid, Indices(0, nEntries),
1087 Values(0, nEntries));
1089 outMatrix->fillComplete();
1091 #ifdef INCLUDE_LENGTHY_OUTPUT
1094 std::cout << inMap->getComm()->getRank() <<
" KDDKDD Rows "
1095 <<
"nGlobal " << inMatrix->getGlobalNumRows() <<
" "
1096 << outMatrix->getGlobalNumRows() <<
"; "
1097 <<
"nLocal " << inMatrix->getLocalNumRows() <<
" "
1098 << outMatrix->getLocalNumRows() << std::endl;
1099 std::cout << inMap->getComm()->getRank() <<
" KDDKDD NNZS "
1100 <<
"nGlobal " << inMatrix->getGlobalNumEntries() <<
" "
1101 << outMatrix->getGlobalNumEntries() <<
"; "
1102 <<
"nLocal " << inMatrix->getLocalNumEntries() <<
" "
1103 << outMatrix->getLocalNumEntries() << std::endl;
1106 Teuchos::Array<zgno_t> in(rowLen), out(rowLen);
1107 Teuchos::Array<zscalar_t> inval(rowLen), outval(rowLen);
1109 for (
size_t i = 0; i < nRows; i++) {
1110 std::cout << inMap->getComm()->getRank() <<
" KDDKDD " << i <<
" nnz(";
1111 inMatrix->getGlobalRowCopy(inMap->getGlobalElement(i), in, inval, nIn);
1112 outMatrix->getGlobalRowCopy(outMap->getGlobalElement(i), out, outval,
1115 std::cout << nIn <<
", " << nOut <<
"): ";
1116 for (
size_t j = 0; j < nIn; j++) {
1117 std::cout <<
"(" << in[j] <<
" " << inval[j] <<
", "
1118 << out[j] <<
" " << outval[j] <<
") ";
1120 std::cout << std::endl;
1123 #endif // INCLUDE_LENGTHY_OUTPUT
1130 void UserInputForTests::readMatrixMarketFile(
1133 bool distributeInput
1136 std::ostringstream
fname;
1137 fname << path <<
"/" << testData <<
".mtx";
1139 if (verbose_ && tcomm_->getRank() == 0)
1140 std::cout <<
"UserInputForTests, Read: " << fname.str() << std::endl;
1145 RCP<tcrsMatrix_t> toMatrix;
1146 RCP<tcrsMatrix_t> fromMatrix;
1149 typedef Tpetra::MatrixMarket::Reader<tcrsMatrix_t> reader_type;
1150 fromMatrix = reader_type::readSparseFile(fname.str(), tcomm_,
1151 true,
false,
false);
1152 #ifdef KDD_NOT_READY_YET
1155 fromMatrix = modifyMatrixGIDs(fromMatrix);
1158 if(!distributeInput)
1160 if (verbose_ && tcomm_->getRank() == 0)
1161 std::cout <<
"Constructing serial distribution of matrix" << std::endl;
1163 RCP<const map_t> fromMap = fromMatrix->getRowMap();
1165 size_t numGlobalCoords = fromMap->getGlobalNumElements();
1166 size_t numLocalCoords = this->tcomm_->getRank() == 0 ? numGlobalCoords : 0;
1167 RCP<const map_t> toMap = rcp(
new map_t(numGlobalCoords,numLocalCoords, 0, tcomm_));
1169 RCP<import_t> importer = rcp(
new import_t(fromMap, toMap));
1171 toMatrix->doImport(*fromMatrix, *importer, Tpetra::INSERT);
1172 toMatrix->fillComplete();
1175 toMatrix = fromMatrix;
1177 }
catch (std::exception &e) {
1178 if (tcomm_->getRank() == 0) {
1179 std::cout <<
"UserInputForTests unable to read matrix market file:"
1180 << fname.str() << std::endl;
1181 std::cout << e.what() << std::endl;
1186 "UserInputForTests unable to read matrix market file");
1189 #ifdef INCLUDE_LENGTHY_OUTPUT
1190 std::cout << tcomm_->getRank() <<
" KDDKDD " << M_->getLocalNumRows()
1191 <<
" " << M_->getGlobalNumRows()
1192 <<
" " << M_->getLocalNumEntries()
1193 <<
" " << M_->getGlobalNumEntries() << std::endl;
1194 #endif // INCLUDE_LENGTHY_OUTPUT
1201 fname << path <<
"/" << testData <<
"_coord.mtx";
1203 size_t coordDim = 0, numGlobalCoords = 0;
1204 size_t msg[2]={0,0};
1205 ArrayRCP<ArrayRCP<zscalar_t> > xyz;
1206 std::ifstream coordFile;
1208 if (tcomm_->getRank() == 0){
1211 std::cout <<
"UserInputForTests, Read: " <<
1212 fname.str() << std::endl;
1216 coordFile.open(fname.str().c_str());
1218 catch (std::exception &){
1229 while (!done && !fail && coordFile.good()){
1230 coordFile.getline(c, 256);
1233 else if (c[0] ==
'%')
1237 std::istringstream s(c);
1238 s >> numGlobalCoords >> coordDim;
1239 if (!s.eof() || numGlobalCoords < 1 || coordDim < 1)
1248 xyz = Teuchos::arcp(
new ArrayRCP<zscalar_t> [coordDim], 0, coordDim);
1250 for (
size_t dim=0; !fail && dim < coordDim; dim++){
1256 xyz[dim] = Teuchos::arcp(tmp, 0, numGlobalCoords);
1258 for (idx=0; !coordFile.eof() && idx < numGlobalCoords; idx++){
1259 coordFile.getline(c, 256);
1260 std::istringstream s(c);
1264 if (idx < numGlobalCoords)
1270 ArrayRCP<zscalar_t> emptyArray;
1271 for (
size_t dim=0; dim < coordDim; dim++)
1272 xyz[dim] = emptyArray;
1285 msg[1] = numGlobalCoords;
1289 Teuchos::broadcast<int, size_t>(*tcomm_, 0, 2, msg);
1292 numGlobalCoords = msg[1];
1298 RCP<const map_t> toMap;
1301 const RCP<const map_t> &mapM = M_->getRowMap();
1305 if (verbose_ && tcomm_->getRank() == 0)
1307 std::cout <<
"Matrix was null. ";
1308 std::cout <<
"Constructing distribution map for coordinate vector."
1312 if(!distributeInput)
1314 if (verbose_ && tcomm_->getRank() == 0)
1315 std::cout <<
"Constructing serial distribution map for coordinates."
1318 size_t numLocalCoords = this->tcomm_->getRank()==0 ? numGlobalCoords : 0;
1319 toMap = rcp(
new map_t(numGlobalCoords,numLocalCoords, base, tcomm_));
1321 toMap = rcp(
new map_t(numGlobalCoords, base, tcomm_));
1329 ArrayRCP<ArrayView<const zscalar_t> > coordLists(coordDim);
1331 if (tcomm_->getRank() == 0){
1333 for (
size_t dim=0; dim < coordDim; dim++)
1334 coordLists[dim] = xyz[dim].view(0, numGlobalCoords);
1338 throw std::bad_alloc();
1340 ArrayRCP<const zgno_t> rowIds = Teuchos::arcp(tmp, 0, numGlobalCoords);
1342 #ifdef KDD_NOT_READY_YET
1351 RCP<const map_t> fromMap = rcp(
new map_t(numGlobalCoords,
1352 rowIds.view(0, numGlobalCoords),
1355 tMVector_t allCoords(fromMap, coordLists.view(0, coordDim), coordDim);
1359 xyz_->doExport(allCoords, exporter, Tpetra::INSERT);
1363 RCP<const map_t> fromMap = rcp(
new map_t(numGlobalCoords,
1364 ArrayView<zgno_t>(), base, tcomm_));
1366 tMVector_t allCoords(fromMap, coordLists.view(0, coordDim), coordDim);
1370 xyz_->doExport(allCoords, exporter, Tpetra::INSERT);
1374 void UserInputForTests::buildCrsMatrix(
int xdim,
int ydim,
int zdim,
1375 string problemType,
bool distributeInput)
1377 #ifdef HAVE_ZOLTAN2_GALERI
1378 Teuchos::CommandLineProcessor tclp;
1379 Galeri::Xpetra::Parameters<zgno_t> params(tclp,
1380 xdim, ydim, zdim, problemType);
1382 RCP<const Tpetra::Map<zlno_t, zgno_t> > map;
1383 if (distributeInput)
1384 map = rcp(
new Tpetra::Map<zlno_t, zgno_t>(params.GetNumGlobalElements(),
1388 size_t nGlobalElements = params.GetNumGlobalElements();
1389 size_t nLocalElements = ((tcomm_->getRank() == 0) ? nGlobalElements : 0);
1390 map = rcp(
new Tpetra::Map<zlno_t, zgno_t>(nGlobalElements, nLocalElements, 0,
1394 if (verbose_ && tcomm_->getRank() == 0){
1396 std::cout <<
"Matrix is " << (distributeInput ?
"" :
"not");
1397 std::cout <<
"distributed." << std::endl;
1399 std::cout <<
"UserInputForTests, Create matrix with " << problemType;
1400 std::cout <<
" (and " << xdim;
1402 std::cout <<
" x " << ydim <<
" x " << zdim;
1404 std::cout <<
" x" << ydim <<
" x 1";
1406 std::cout <<
"x 1 x 1";
1408 std::cout <<
" mesh)" << std::endl;
1414 RCP<Galeri::Xpetra::Problem<Tpetra::Map<zlno_t, zgno_t>, Tpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t>, Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t> > > Pr =
1415 Galeri::Xpetra::BuildProblem<zscalar_t, zlno_t, zgno_t, Tpetra::Map<zlno_t, zgno_t>, Tpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t>, Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t> >
1416 (params.GetMatrixType(), map, params.GetParameterList());
1417 M_ = Pr->BuildMatrix();
1419 catch (std::exception &) {
1423 "UserInputForTests Galeri::Xpetra::BuildProblem failed");
1429 if (verbose_ && tcomm_->getRank() == 0)
1431 "UserInputForTests, Implied matrix row coordinates computed" <<
1434 ArrayView<const zgno_t> gids = map->getLocalElementList();
1437 size_t pos = problemType.find(
"2D");
1438 if (pos != string::npos)
1440 else if (problemType ==
string(
"Laplace1D") ||
1441 problemType ==
string(
"Identity"))
1447 for (
int i=0; i < dim; i++){
1450 throw(std::bad_alloc());
1451 coordinates[i] = Teuchos::arcp(c, 0, count,
true);
1458 zgno_t xySize = xdim * ydim;
1459 for (
zlno_t i=0; i < count; i++){
1460 zgno_t iz = gids[i] / xySize;
1461 zgno_t xy = gids[i] - iz*xySize;
1470 for (
zlno_t i=0; i < count; i++){
1477 for (
zlno_t i=0; i < count; i++)
1482 Array<ArrayView<const zscalar_t> > coordView(dim);
1484 for (
int i=0; i < dim; i++)
1487 xyz_ = rcp(
new tMVector_t(map, coordView.view(0, dim), dim));
1489 throw std::runtime_error(
"Galeri input requested but Trilinos is "
1490 "not built with Galeri.");
1494 void UserInputForTests::readZoltanTestData(
string path,
string testData,
1495 bool distributeInput)
1497 int rank = tcomm_->getRank();
1498 FILE *graphFile = NULL;
1499 FILE *coordFile = NULL;
1500 FILE *assignFile = NULL;
1503 for (
int i = 0; i < CHACO_LINE_LENGTH; i++) chaco_line[i] =
'\0';
1507 std::ostringstream chGraphFileName;
1508 chGraphFileName << path <<
"/" << testData <<
".graph";
1511 std::ostringstream chCoordFileName;
1512 chCoordFileName << path <<
"/" << testData <<
".coords";
1515 std::ostringstream chAssignFileName;
1516 chAssignFileName << path <<
"/" << testData <<
".assign";
1519 graphFile = fopen(chGraphFileName.str().c_str(),
"r");
1523 chGraphFileName.str(
"");
1524 chCoordFileName.str(
"");
1526 chGraphFileName << path <<
"/ch_" << testData <<
"/" << testData <<
".graph";
1527 chCoordFileName << path <<
"/ch_" << testData <<
"/" << testData <<
".coords";
1528 chAssignFileName << path <<
"/ch_" << testData <<
"/" << testData <<
".assign";
1531 graphFile = fopen(chGraphFileName.str().c_str(),
"r");
1534 memset(fileInfo, 0,
sizeof(
int) * 3);
1537 if (verbose_ && tcomm_->getRank() == 0)
1538 std::cout <<
"UserInputForTests, open " <<
1539 chGraphFileName.str () << std::endl;
1541 coordFile = fopen(chCoordFileName.str().c_str(),
"r");
1544 if (verbose_ && tcomm_->getRank() == 0)
1545 std::cout <<
"UserInputForTests, open " <<
1546 chCoordFileName.str () << std::endl;
1549 assignFile = fopen(chAssignFileName.str().c_str(),
"r");
1552 if (verbose_ && tcomm_->getRank() == 0)
1553 std::cout <<
"UserInputForTests, open " <<
1554 chAssignFileName.str () << std::endl;
1557 if (verbose_ && tcomm_->getRank() == 0){
1558 std::cout <<
"UserInputForTests, unable to open file: ";
1559 std::cout << chGraphFileName.str() << std::endl;
1565 Teuchos::broadcast<int, int>(*tcomm_, 0, 3, fileInfo);
1567 bool haveGraph = (fileInfo[0] == 1);
1568 bool haveCoords = (fileInfo[1] == 1);
1569 bool haveAssign = (fileInfo[2] == 1);
1574 getUIChacoGraph(graphFile, haveAssign, assignFile,
1575 testData, distributeInput);
1582 getUIChacoCoords(coordFile, testData);
1591 void UserInputForTests::getUIChacoGraph(FILE *fptr,
bool haveAssign,
1592 FILE *assignFile,
string fname,
1593 bool distributeInput)
1595 int rank = tcomm_->getRank();
1597 int nvtxs=0, nedges=0;
1598 int nVwgts=0, nEwgts=0;
1599 int *start = NULL, *adj = NULL;
1600 float *ewgts = NULL, *vwgts = NULL;
1601 size_t *nzPerRow = NULL;
1602 size_t maxRowLen = 0;
1604 ArrayRCP<const size_t> rowSizes;
1606 bool haveEdges =
true;
1610 memset(graphCounts, 0, 5*
sizeof(
int));
1613 fail = chaco_input_graph(fptr, fname.c_str(), &start, &adj,
1614 &nvtxs, &nVwgts, &vwgts, &nEwgts, &ewgts);
1625 std::cout <<
"UserInputForTests, " << nvtxs <<
" vertices,";
1627 std::cout << start[nvtxs] <<
" edges,";
1629 std::cout <<
"no edges,";
1630 std::cout << nVwgts <<
" vertex weights, ";
1631 std::cout << nEwgts <<
" edge weights" << std::endl;
1638 Teuchos::broadcast<int, int>(*tcomm_, 0, 5, graphCounts);
1639 throw std::runtime_error(
"Unable to read chaco file");
1643 nedges = start[nvtxs];
1645 nzPerRow =
new size_t [nvtxs];
1647 throw std::bad_alloc();
1648 rowSizes = arcp(nzPerRow, 0, nvtxs,
true);
1651 for (
int i=0; i < nvtxs; i++){
1652 nzPerRow[i] = start[i+1] - start[i];
1653 if (nzPerRow[i] > maxRowLen)
1654 maxRowLen = nzPerRow[i];
1658 memset(nzPerRow, 0,
sizeof(
size_t) * nvtxs);
1665 for (
int i=0; i < nedges; i++)
1669 graphCounts[0] = nvtxs;
1670 graphCounts[1] = nedges;
1671 graphCounts[2] = nVwgts;
1672 graphCounts[3] = nEwgts;
1673 graphCounts[4] = (int)maxRowLen;
1676 Teuchos::broadcast<int, int>(*tcomm_, 0, 5, graphCounts);
1678 if (graphCounts[0] == 0)
1679 throw std::runtime_error(
"Unable to read chaco file");
1681 haveEdges = (graphCounts[1] > 0);
1683 RCP<tcrsMatrix_t> fromMatrix;
1684 RCP<const map_t> fromMap;
1688 fromMap = rcp(
new map_t(nvtxs, nvtxs, base, tcomm_));
1696 if (nedges && !edgeIds)
1697 throw std::bad_alloc();
1698 for (
int i=0; i < nedges; i++)
1699 edgeIds[i] = adj[i];
1704 zgno_t *nextId = edgeIds;
1705 Array<zscalar_t> values(nedges, 1.0);
1706 if (nedges > 0 && nEwgts > 0) {
1707 for (
int i=0; i < nedges; i++)
1708 values[i] = ewgts[i];
1713 for (
int i=0; i < nvtxs; i++){
1714 if (nzPerRow[i] > 0){
1715 ArrayView<const zgno_t> rowNz(nextId, nzPerRow[i]);
1716 fromMatrix->insertGlobalValues(i, rowNz, values.view(start[i], start[i+1] - start[i]));
1717 nextId += nzPerRow[i];
1725 fromMatrix->fillComplete();
1728 nvtxs = graphCounts[0];
1729 nedges = graphCounts[1];
1730 nVwgts = graphCounts[2];
1731 nEwgts = graphCounts[3];
1732 maxRowLen = graphCounts[4];
1736 fromMap = rcp(
new map_t(nvtxs, 0, base, tcomm_));
1741 fromMatrix->fillComplete();
1746 size_t sz = fromMatrix->getLocalMaxNumRowEntries();
1747 Teuchos::Array<zgno_t> indices(sz);
1748 Teuchos::Array<zscalar_t> values(sz);
1749 for (
size_t i = 0; i < fromMatrix->getLocalNumRows(); i++) {
1750 zgno_t gid = fromMatrix->getRowMap()->getGlobalElement(i);
1752 fromMatrix->getGlobalRowCopy(gid, indices(), values(), num);
1753 std::cout <<
"ROW " << gid <<
": ";
1754 for (
size_t j = 0; j < num; j++)
1755 std::cout << indices[j] <<
" ";
1756 std::cout << std::endl;
1761 RCP<const map_t> toMap;
1762 RCP<tcrsMatrix_t> toMatrix;
1763 RCP<import_t> importer;
1765 if (distributeInput) {
1768 short *assignments =
new short[nvtxs];
1770 fail = chaco_input_assign(assignFile, fname.c_str(), nvtxs, assignments);
1773 Teuchos::broadcast<int, short>(*tcomm_, 0, nvtxs, assignments);
1776 Teuchos::Array<zgno_t> mine;
1777 for (
int i = 0; i < nvtxs; i++) {
1778 if (assignments[i] == rank)
1783 Teuchos::OrdinalTraits<Tpetra::global_size_t>::invalid();
1784 toMap = rcp(
new map_t(dummy, mine(), base, tcomm_));
1785 delete [] assignments;
1789 toMap = rcp(
new map_t(nvtxs, base, tcomm_));
1794 importer = rcp(
new import_t(fromMap, toMap));
1795 toMatrix->doImport(*fromMatrix, *importer, Tpetra::INSERT);
1796 toMatrix->fillComplete();
1800 toMatrix = fromMatrix;
1807 typedef ArrayRCP<const ArrayView<const zscalar_t> > arrayArray_t;
1811 ArrayRCP<zscalar_t> weightBuf;
1812 ArrayView<const zscalar_t> *wgts =
new ArrayView<const zscalar_t> [nVwgts];
1815 size_t len = nVwgts * nvtxs;
1817 if (!buf)
throw std::bad_alloc();
1818 weightBuf = arcp(buf, 0, len,
true);
1820 for (
int widx=0; widx < nVwgts; widx++){
1821 wgts[widx] = ArrayView<const zscalar_t>(buf, nvtxs);
1822 float *vw = vwgts + widx;
1823 for (
int i=0; i < nvtxs; i++, vw += nVwgts)
1832 arrayArray_t vweights = arcp(wgts, 0, nVwgts,
true);
1834 RCP<tMVector_t> fromVertexWeights =
1835 rcp(
new tMVector_t(fromMap, vweights.view(0, nVwgts), nVwgts));
1837 RCP<tMVector_t> toVertexWeights;
1838 if (distributeInput) {
1839 toVertexWeights = rcp(
new tMVector_t(toMap, nVwgts));
1840 toVertexWeights->doImport(*fromVertexWeights, *importer, Tpetra::INSERT);
1843 toVertexWeights = fromVertexWeights;
1845 vtxWeights_ = toVertexWeights;
1850 if (haveEdges && nEwgts > 0){
1900 toMap = rcp(
new map_t(nedges, M_->getLocalNumEntries(), base, tcomm_));
1901 edgWeights_ = rcp(
new tMVector_t(toMap, nEwgts));
1903 size_t maxSize = M_->getLocalMaxNumRowEntries();
1904 tcrsMatrix_t::nonconst_local_inds_host_view_type colind(
"colind", maxSize);
1905 tcrsMatrix_t::nonconst_values_host_view_type vals(
"values", maxSize);
1908 for (
size_t i = 0, idx = 0; i < M_->getLocalNumRows(); i++) {
1909 M_->getLocalRowCopy(i, colind, vals, nEntries);
1910 for (
size_t j = 0; j < nEntries; j++) {
1911 edgWeights_->replaceLocalValue(idx, 0, vals[j]);
1923 void UserInputForTests::getUIChacoCoords(FILE *fptr,
string fname)
1925 int rank = tcomm_->getRank();
1927 double *x=NULL, *y=NULL, *z=NULL;
1930 size_t globalNumVtx = M_->getGlobalNumRows();
1937 fail = chaco_input_geom(fptr, fname.c_str(), (int)globalNumVtx,
1944 std::cout <<
"UserInputForTests, read " << globalNumVtx;
1945 std::cout <<
" " << ndim <<
"-dimensional coordinates." << std::endl;
1949 Teuchos::broadcast<int, int>(*tcomm_, 0, 1, &ndim);
1952 throw std::runtime_error(
"Can't read coordinate file");
1954 ArrayRCP<ArrayRCP<const zscalar_t> > coords(ndim);
1959 for (
int dim=0; dim < ndim; dim++){
1962 throw std::bad_alloc();
1963 coords[dim] = arcp<const zscalar_t>(v, 0, globalNumVtx,
true);
1964 double *val = (dim==0 ? x : (dim==1 ? y : z));
1965 for (
size_t i=0; i < globalNumVtx; i++)
1971 len =
static_cast<zlno_t>(globalNumVtx);
1974 RCP<const map_t> fromMap = rcp(
new map_t(globalNumVtx, len, 0, tcomm_));
1975 RCP<const map_t> toMap = M_->getRowMap();
1976 RCP<import_t> importer = rcp(
new import_t(fromMap, toMap));
1978 Array<ArrayView<const zscalar_t> > coordData;
1979 for (
int dim=0; dim < ndim; dim++)
1980 coordData.push_back(coords[dim].view(0, len));
1982 RCP<tMVector_t> fromCoords =
1983 rcp(
new tMVector_t(fromMap, coordData.view(0, ndim), ndim));
1985 RCP<tMVector_t> toCoords = rcp(
new tMVector_t(toMap, ndim));
1987 toCoords->doImport(*fromCoords, *importer, Tpetra::INSERT);
1996 double UserInputForTests::chaco_read_val(
2012 if (chaco_offset == 0 || chaco_offset >= chaco_break_pnt) {
2013 if (chaco_offset >= chaco_break_pnt) {
2014 length_left = CHACO_LINE_LENGTH - chaco_save_pnt - 1;
2016 ptr = &chaco_line[chaco_save_pnt];
2017 for (i=length_left; i; i--) *ptr2++ = *ptr++;
2018 length = chaco_save_pnt + 1;
2021 length = CHACO_LINE_LENGTH;
2026 ptr2 = fgets(&chaco_line[length_left], length, infile);
2028 if (ptr2 == (
char *) NULL) {
2030 return((
double) 0.0);
2033 if ((chaco_line[CHACO_LINE_LENGTH - 2] !=
'\n') && (chaco_line[CHACO_LINE_LENGTH - 2] !=
'\f')
2034 && (strlen(chaco_line) == CHACO_LINE_LENGTH - 1)){
2036 chaco_break_pnt = CHACO_LINE_LENGTH - 1;
2037 chaco_save_pnt = chaco_break_pnt;
2042 if (chaco_line[chaco_break_pnt] !=
'\0') {
2043 if (isspace((
int)(chaco_line[chaco_break_pnt]))) {
2045 chaco_save_pnt = chaco_break_pnt + 1;
2049 else if (white_seen) {
2056 chaco_break_pnt = CHACO_LINE_LENGTH;
2062 while (isspace((
int)(chaco_line[chaco_offset])) && chaco_offset < CHACO_LINE_LENGTH) chaco_offset++;
2063 if (chaco_line[chaco_offset] ==
'%' || chaco_line[chaco_offset] ==
'#') {
2065 if (chaco_break_pnt < CHACO_LINE_LENGTH) {
2066 chaco_flush_line(infile);
2068 return((
double) 0.0);
2071 ptr = &(chaco_line[chaco_offset]);
2072 val = strtod(ptr, &ptr2);
2077 return((
double) 0.0);
2080 chaco_offset = (int) (ptr2 - chaco_line) /
sizeof(char);
2087 int UserInputForTests::chaco_read_int(
2103 if (chaco_offset == 0 || chaco_offset >= chaco_break_pnt) {
2104 if (chaco_offset >= chaco_break_pnt) {
2105 length_left = CHACO_LINE_LENGTH - chaco_save_pnt - 1;
2107 ptr = &chaco_line[chaco_save_pnt];
2108 for (i=length_left; i; i--) *ptr2++ = *ptr++;
2109 length = chaco_save_pnt + 1;
2112 length = CHACO_LINE_LENGTH;
2117 ptr2 = fgets(&chaco_line[length_left], length, infile);
2119 if (ptr2 == (
char *) NULL) {
2124 if ((chaco_line[CHACO_LINE_LENGTH - 2] !=
'\n') && (chaco_line[CHACO_LINE_LENGTH - 2] !=
'\f')
2125 && (strlen(chaco_line) == CHACO_LINE_LENGTH - 1)){
2127 chaco_break_pnt = CHACO_LINE_LENGTH - 1;
2128 chaco_save_pnt = chaco_break_pnt;
2133 if (chaco_line[chaco_break_pnt] !=
'\0') {
2134 if (isspace((
int)(chaco_line[chaco_break_pnt]))) {
2136 chaco_save_pnt = chaco_break_pnt + 1;
2140 else if (white_seen) {
2147 chaco_break_pnt = CHACO_LINE_LENGTH;
2153 while (isspace((
int)(chaco_line[chaco_offset])) && chaco_offset < CHACO_LINE_LENGTH) chaco_offset++;
2154 if (chaco_line[chaco_offset] ==
'%' || chaco_line[chaco_offset] ==
'#') {
2156 if (chaco_break_pnt < CHACO_LINE_LENGTH) {
2157 chaco_flush_line(infile);
2162 ptr = &(chaco_line[chaco_offset]);
2163 val = (int) strtol(ptr, &ptr2, 10);
2171 chaco_offset = (int) (ptr2 - chaco_line) /
sizeof(char);
2177 void UserInputForTests::chaco_flush_line(
2184 while (c !=
'\n' && c !=
'\f')
2188 int UserInputForTests::chaco_input_graph(
2238 while (end_flag == 1) {
2239 *nvtxs = chaco_read_int(fin, &end_flag);
2243 printf(
"ERROR in graph file `%s':", inname);
2244 printf(
" Invalid number of vertices (%d).\n", *nvtxs);
2249 narcs = chaco_read_int(fin, &end_flag);
2251 printf(
"ERROR in graph file `%s':", inname);
2252 printf(
" Invalid number of expected edges (%d).\n", narcs);
2259 option = chaco_read_int(fin, &end_flag);
2261 using_ewgts = option - 10 * (option / 10);
2263 using_vwgts = option - 10 * (option / 10);
2265 vtxnums = option - 10 * (option / 10);
2268 (*nVwgts) = using_vwgts;
2269 (*nEwgts) = using_ewgts;
2272 if (!end_flag && using_vwgts==1){
2273 j = chaco_read_int(fin, &end_flag);
2274 if (!end_flag) (*nVwgts) = j;
2276 if (!end_flag && using_ewgts==1){
2277 j = chaco_read_int(fin, &end_flag);
2278 if (!end_flag) (*nEwgts) = j;
2283 j = chaco_read_int(fin, &end_flag);
2286 *start = (
int *) malloc((
unsigned) (*nvtxs + 1) *
sizeof(
int));
2288 *adjacency = (
int *) malloc((
unsigned) (2 * narcs + 1) *
sizeof(
int));
2293 *vweights = (
float *) malloc((
unsigned) (*nvtxs) * (*nVwgts) *
sizeof(float));
2298 *eweights = (
float *)
2299 malloc((
unsigned) (2 * narcs + 1) * (*nEwgts) *
sizeof(float));
2303 adjptr = *adjacency;
2312 while ((using_vwgts || vtxnums || narcs) && end_flag != -1) {
2317 j = chaco_read_int(fin, &end_flag);
2319 if (vertex == *nvtxs)
2321 printf(
"ERROR in graph file `%s':", inname);
2322 printf(
" no vertex number in line %d.\n", line_num);
2326 if (j != vertex && j != vertex + 1) {
2327 printf(
"ERROR in graph file `%s':", inname);
2328 printf(
" out-of-order vertex number in line %d.\n", line_num);
2342 if (vertex > *nvtxs)
2346 if (using_vwgts && new_vertex) {
2347 for (j=0; j<(*nVwgts); j++){
2348 weight = chaco_read_val(fin, &end_flag);
2350 printf(
"ERROR in graph file `%s':", inname);
2351 printf(
" not enough weights for vertex %d.\n", vertex);
2355 (*vweights)[(vertex-1)*(*nVwgts)+j] = weight;
2362 neighbor = chaco_read_int(fin, &end_flag);
2368 for (j=0; j<(*nEwgts); j++){
2369 eweight = chaco_read_val(fin, &end_flag);
2372 printf(
"ERROR in graph file `%s':", inname);
2373 printf(
" not enough weights for edge (%d,%d).\n", vertex, neighbor);
2386 if (++nedges > 2*narcs) {
2387 printf(
"ERROR in graph file `%s':", inname);
2388 printf(
" at least %d adjacencies entered, but nedges = %d\n",
2393 *adjptr++ = neighbor;
2398 neighbor = chaco_read_int(fin, &end_flag);
2402 (*start)[vertex] = sum_edges;
2407 while (!flag && end_flag != -1) {
2408 chaco_read_int(fin, &end_flag);
2413 (*start)[*nvtxs] = sum_edges;
2421 if (*adjacency != NULL)
2424 if (*eweights != NULL)
2433 if (*adjacency != NULL)
2435 if (*vweights != NULL)
2437 if (*eweights != NULL)
2445 return (error_flag);
2449 int UserInputForTests::chaco_input_geom(
2451 const char *geomname,
2459 double xc, yc, zc =0;
2467 *x = *y = *z = NULL;
2470 while (end_flag == 1) {
2471 xc = chaco_read_val(fingeom, &end_flag);
2475 if (end_flag == -1) {
2476 printf(
"No values found in geometry file `%s'\n", geomname);
2482 yc = chaco_read_val(fingeom, &end_flag);
2483 if (end_flag == 0) {
2485 zc = chaco_read_val(fingeom, &end_flag);
2486 if (end_flag == 0) {
2488 chaco_read_val(fingeom, &end_flag);
2490 printf(
"Too many values on input line of geometry file `%s'\n",
2493 printf(
" Maximum dimensionality is 3\n");
2502 *x = (
double *) malloc((
unsigned) nvtxs *
sizeof(double));
2505 *y = (
double *) malloc((
unsigned) nvtxs *
sizeof(double));
2509 *z = (
double *) malloc((
unsigned) nvtxs *
sizeof(double));
2513 for (nread = 1; nread < nvtxs; nread++) {
2516 i = fscanf(fingeom,
"%lf", &((*x)[nread]));
2518 else if (ndims == 2) {
2519 i = fscanf(fingeom,
"%lf%lf", &((*x)[nread]), &((*y)[nread]));
2521 else if (ndims == 3) {
2522 i = fscanf(fingeom,
"%lf%lf%lf", &((*x)[nread]), &((*y)[nread]),
2527 printf(
"Too few lines of values in geometry file; nvtxs=%d, but only %d read\n",
2532 else if (i != ndims) {
2533 printf(
"Wrong number of values in line %d of geometry file `%s'\n",
2534 line_num, geomname);
2543 while (!flag && end_flag != -1) {
2544 chaco_read_val(fingeom, &end_flag);
2556 int UserInputForTests::chaco_input_assign(
2558 const char *inassignname,
2569 while (end_flag == 1) {
2570 assignment[0] = chaco_read_int(finassign, &end_flag);
2573 if (assignment[0] < 0) {
2574 printf(
"ERROR: Entry %d in assignment file `%s' less than zero (%d)\n",
2575 1, inassignname, assignment[0]);
2580 if (end_flag == -1) {
2581 printf(
"ERROR: No values found in assignment file `%s'\n", inassignname);
2587 int np = this->tcomm_->getSize();
2588 if (assignment[0] >= np) flag = assignment[0];
2589 srand(this->tcomm_->getRank());
2590 for (i = 1; i < nvtxs; i++) {
2591 j = fscanf(finassign,
"%hd", &(assignment[i]));
2593 printf(
"ERROR: Too few values in assignment file `%s'.\n", inassignname);
2597 if (assignment[i] < 0) {
2598 printf(
"ERROR: Entry %d in assignment file `%s' less than zero (%d)\n",
2599 i+1, inassignname, assignment[i]);
2603 if (assignment[i] >= np) {
2605 if (assignment[i] > flag)
2606 flag = assignment[i];
2607 assignment[i] = rand() % np;
2613 printf(
"WARNING: Possible error in assignment file `%s'\n",
2615 printf(
" Max assignment set (%d) greater than "
2616 "max processor rank (%d)\n", flag, np-1);
2617 printf(
" Some vertices given random initial assignments\n");
2623 while (!flag && end_flag != -1) {
2624 chaco_read_int(finassign, &end_flag);
2629 printf(
"WARNING: Possible error in assignment file `%s'\n", inassignname);
2630 printf(
" Numerical data found after expected end of file\n");
2638 void UserInputForTests::readPamgenMeshFile(
string path,
string testData)
2640 #ifdef HAVE_ZOLTAN2_PAMGEN
2641 int rank = this->tcomm_->getRank();
2642 if (verbose_ && tcomm_->getRank() == 0)
2643 std::cout <<
"UserInputForTestsBD::readPamgenFile, Read: " << testData << std::endl;
2650 std::ostringstream meshFileName;
2651 meshFileName << path <<
"/" << testData <<
".pmgen";
2654 file.open(meshFileName.str(), std::ios::in);
2658 if(verbose_ && tcomm_->getRank() == 0)
2660 std::cout <<
"Unable to open pamgen mesh: ";
2661 std::cout << meshFileName.str();
2662 std::cout <<
"\nPlease check file path and name." << std::endl;
2668 file.seekg (0,file.end);
2675 while(std::getline(file,line))
2677 if( line.find(
"nz") != std::string::npos ||
2678 line.find(
"nphi") != std::string::npos)
2686 file.seekg(0, std::ios::beg);
2691 this->tcomm_->broadcast(0,
sizeof(
int), (
char *)&dimension);
2692 this->tcomm_->broadcast(0,
sizeof(
size_t),(
char *)&len);
2693 this->tcomm_->barrier();
2696 if(verbose_ && tcomm_->getRank() == 0)
2697 std::cout <<
"Pamgen Mesh file size == 0, exiting UserInputForTests early." << std::endl;
2701 char * file_data =
new char[len+1];
2702 file_data[len] =
'\0';
2704 file.read(file_data,len);
2708 this->tcomm_->broadcast(0,(
int)len+1,file_data);
2709 this->tcomm_->barrier();
2713 this->pamgen_mesh = rcp(
new PamgenMesh);
2714 this->havePamgenMesh =
true;
2715 pamgen_mesh->createMesh(file_data,dimension,this->tcomm_);
2718 pamgen_mesh->storeMesh();
2719 this->tcomm_->barrier();
2722 this->setPamgenCoordinateMV();
2725 this->setPamgenAdjacencyGraph();
2727 this->tcomm_->barrier();
2728 if(rank == 0) file.close();
2729 delete [] file_data;
2731 throw std::runtime_error(
"Pamgen requested but Trilinos "
2732 "not built with Pamgen");
2736 #ifdef HAVE_ZOLTAN2_PAMGEN
2737 void UserInputForTests::setPamgenCoordinateMV()
2739 int dimension = pamgen_mesh->num_dim;
2743 zgno_t numelements = pamgen_mesh->num_elem;
2744 zgno_t numGlobalElements = pamgen_mesh->num_elems_global;
2747 for(
int i = 0; i < dimension; ++i){
2748 elem_coords[i] =
new zscalar_t[numelements];
2749 double *tmp = &pamgen_mesh->element_coord[i*numelements];
2750 for (
int j = 0; j < numelements; j++) elem_coords[i][j] = tmp[j];
2754 RCP<Tpetra::Map<zlno_t, zgno_t, znode_t> > mp;
2758 Array<zgno_t> elementList(numelements);
2759 for (Array<zgno_t>::size_type k = 0; k < numelements; ++k) {
2760 elementList[k] = pamgen_mesh->element_order_map[k];
2763 mp = rcp (
new map_t (numGlobalElements, elementList, 0, this->tcomm_));
2767 Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > coordView(dimension);
2768 for (
int i = 0; i < dimension; i++){
2769 if(numelements > 0){
2770 Teuchos::ArrayView<const zscalar_t> a(elem_coords[i], numelements);
2774 Teuchos::ArrayView<const zscalar_t> a;
2780 xyz_ = RCP<tMVector_t>(
new
2785 void UserInputForTests::setPamgenAdjacencyGraph()
2795 size_t local_nodes = (size_t)this->pamgen_mesh->num_nodes;
2796 size_t local_els = (
size_t)this->pamgen_mesh->num_elem;
2798 size_t global_els = (size_t)this->pamgen_mesh->num_elems_global;
2799 size_t global_nodes = (
size_t)this->pamgen_mesh->num_nodes_global;
2803 RCP<const map_t> rowMap = rcp(
new map_t(global_els,0,this->tcomm_));
2804 RCP<const map_t> rangeMap = rowMap;
2807 RCP<const map_t> domainMap = rcp(
new map_t(global_nodes,0,this->tcomm_));
2810 int blks = this->pamgen_mesh->num_elem_blk;
2811 int max_nodes_per_el = 0;
2812 for(
int i = 0; i < blks; i++)
2813 if (this->pamgen_mesh->nodes_per_element[i] > max_nodes_per_el)
2814 max_nodes_per_el = this->pamgen_mesh->nodes_per_element[i];
2817 Teuchos::RCP<tcrsMatrix_t> C = rcp(
new tcrsMatrix_t(rowMap,max_nodes_per_el));
2819 Array<zgno_t> g_el_ids(local_els);
2820 for (
size_t k = 0; k < local_els; ++k) {
2821 g_el_ids[k] = pamgen_mesh->global_element_numbers[k]-1;
2824 Array<zgno_t> g_node_ids(local_nodes);
2825 for (
size_t k = 0; k < local_nodes; ++k) {
2826 g_node_ids[k] = pamgen_mesh->global_node_numbers[k]-1;
2833 for(
int i = 0; i < blks; i++)
2835 int el_per_block = this->pamgen_mesh->elements[i];
2836 int nodes_per_el = this->pamgen_mesh->nodes_per_element[i];
2837 int * connect = this->pamgen_mesh->elmt_node_linkage[i];
2839 for(
int j = 0; j < el_per_block; j++)
2841 const zgno_t gid =
static_cast<zgno_t>(g_el_ids[el_no]);
2842 for(
int k = 0; k < nodes_per_el; k++)
2844 int g_node_i = g_node_ids[connect[j*nodes_per_el+k]-1];
2845 C->insertGlobalValues(gid,
2846 Teuchos::tuple<zgno_t>(g_node_i),
2847 Teuchos::tuple<zscalar_t>(one));
2853 C->fillComplete(domainMap, rangeMap);
2859 Tpetra::MatrixMatrix::Multiply(*C,
false, *C,
true, *A);
2864 this->M_ = rcp(
new tcrsMatrix_t(rowMap, A->getGlobalMaxNumRowEntries()));
2867 Teuchos::ArrayView<const zgno_t> rowMapElementList =
2868 rowMap->getLocalElementList();
2869 for (Teuchos_Ordinal ii = 0; ii < rowMapElementList.size(); ii++)
2871 zgno_t gid = rowMapElementList[ii];
2872 size_t numEntriesInRow = A->getNumEntriesInGlobalRow (gid);
2873 typename tcrsMatrix_t::nonconst_global_inds_host_view_type rowinds(
"Indices", numEntriesInRow);
2874 typename tcrsMatrix_t::nonconst_values_host_view_type rowvals(
"Values", numEntriesInRow);
2877 Array<zscalar_t> mod_rowvals;
2878 Array<zgno_t> mod_rowinds;
2879 A->getGlobalRowCopy (gid, rowinds, rowvals, numEntriesInRow);
2880 for (
size_t i = 0; i < numEntriesInRow; i++) {
2883 if (rowvals[i] >= 1)
2885 mod_rowvals.push_back(one);
2886 mod_rowinds.push_back(rowinds[i]);
2889 this->M_->insertGlobalValues(gid, mod_rowinds, mod_rowvals);
2892 this->M_->fillComplete();
keep typedefs that commonly appear in many places localized
#define TEST_FAIL_AND_THROW(comm, ok, s)
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
string trim_right_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
void readGeoGenParams(string paramFileName, Teuchos::ParameterList &geoparams, const RCP< const Teuchos::Comm< int > > &comm)
Traits of Xpetra classes, including migration method.
int buildCrsMatrix(int xdim, int ydim, int zdim, std::string problemType, const Teuchos::RCP< const Teuchos::Comm< int > > &comm, Teuchos::RCP< Tpetra::CrsMatrix< scalar_t, lno_t, gno_t, nod_t >> &M_)
common code used by tests
static RCP< User > convertToXpetra(const RCP< User > &a)
Convert the object to its Xpetra wrapped version.
string trim_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
GeometricGen::GeometricGenerator< zscalar_t, zlno_t, zgno_t, znode_t > geometricgen_t
Xpetra::CrsGraph< zlno_t, zgno_t, znode_t > xcrsGraph_t
Tpetra::Vector< zscalar_t, zlno_t, zgno_t, znode_t > tVector_t
Tpetra::Map::local_ordinal_type zlno_t
static const std::string fail
Tpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > tcrsMatrix_t
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > tcrsGraph_t
string trim_left_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
Tpetra::global_size_t global_size_t
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tMVector_t
Tpetra::Map::global_ordinal_type zgno_t