50 #ifndef USERINPUTFORTESTS 51 #define USERINPUTFORTESTS 56 #include <Tpetra_MultiVector.hpp> 57 #include <Tpetra_CrsMatrix.hpp> 58 #include <Tpetra_Map.hpp> 59 #include <Xpetra_Vector.hpp> 60 #include <Xpetra_CrsMatrix.hpp> 61 #include <Xpetra_CrsGraph.hpp> 63 #include <MatrixMarket_Tpetra.hpp> 65 #ifdef HAVE_ZOLTAN2_GALERI 66 #include <Galeri_XpetraProblemFactory.hpp> 67 #include <Galeri_XpetraParameters.hpp> 70 #include <Kokkos_DefaultNode.hpp> 76 #include <TpetraExt_MatrixMatrix_def.hpp> 83 using Teuchos::ArrayRCP;
84 using Teuchos::ArrayView;
89 using Teuchos::rcp_const_cast;
90 using Teuchos::ParameterList;
129 typedef Tpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t>
tcrsMatrix_t;
131 typedef Tpetra::Vector<zscalar_t, zlno_t, zgno_t, znode_t>
tVector_t;
132 typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t>
tMVector_t;
134 typedef Xpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t>
xcrsMatrix_t;
136 typedef Xpetra::Vector<zscalar_t, zlno_t, zgno_t, znode_t>
xVector_t;
137 typedef Xpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t>
xMVector_t;
139 typedef Tpetra::Map<zlno_t, zgno_t, znode_t>
map_t;
140 typedef Tpetra::Export<zlno_t, zgno_t, znode_t>
export_t;
141 typedef Tpetra::Import<zlno_t, zgno_t, znode_t>
import_t;
167 const RCP<
const Comm<int> > &c,
bool debugInfo=
false,
168 bool distributeInput=
true);
187 const RCP<
const Comm<int> > &c,
bool debugInfo=
false,
188 bool distributeInput=
true);
206 const RCP<
const Comm<int> > &c);
210 static void getUIRandomData(
unsigned int seed,
zlno_t length,
213 RCP<tMVector_t> getUICoordinates();
215 RCP<tMVector_t> getUIWeights();
217 RCP<tMVector_t> getUIEdgeWeights();
219 RCP<tcrsMatrix_t> getUITpetraCrsMatrix();
221 RCP<tcrsGraph_t> getUITpetraCrsGraph();
223 RCP<tVector_t> getUITpetraVector();
225 RCP<tMVector_t> getUITpetraMultiVector(
int nvec);
227 RCP<xcrsMatrix_t> getUIXpetraCrsMatrix();
229 RCP<xcrsGraph_t> getUIXpetraCrsGraph();
231 RCP<xVector_t> getUIXpetraVector();
233 RCP<xMVector_t> getUIXpetraMultiVector(
int nvec);
235 #ifdef HAVE_ZOLTAN2_PAMGEN 236 PamgenMesh * getPamGenMesh(){
return this->pamgen_mesh.operator->();}
239 #ifdef HAVE_EPETRA_DATA_TYPES 240 RCP<Epetra_CrsGraph> getUIEpetraCrsGraph();
242 RCP<Epetra_CrsMatrix> getUIEpetraCrsMatrix();
244 RCP<Epetra_Vector> getUIEpetraVector();
246 RCP<Epetra_MultiVector> getUIEpetraMultiVector(
int nvec);
250 bool hasInputDataType(
const string &input_type);
252 bool hasUICoordinates();
256 bool hasUIEdgeWeights();
258 bool hasUITpetraCrsMatrix();
260 bool hasUITpetraCrsGraph();
262 bool hasUITpetraVector();
264 bool hasUITpetraMultiVector();
266 bool hasUIXpetraCrsMatrix();
268 bool hasUIXpetraCrsGraph();
270 bool hasUIXpetraVector();
272 bool hasUIXpetraMultiVector();
274 bool hasPamgenMesh();
275 #ifdef HAVE_EPETRA_DATA_TYPES 276 bool hasUIEpetraCrsGraph();
278 bool hasUIEpetraCrsMatrix();
280 bool hasUIEpetraVector();
282 bool hasUIEpetraMultiVector();
290 const RCP<const Comm<int> > tcomm_;
293 #ifdef HAVE_ZOLTAN2_PAMGEN 294 RCP<PamgenMesh> pamgen_mesh;
297 RCP<tcrsMatrix_t> M_;
298 RCP<xcrsMatrix_t> xM_;
300 RCP<tMVector_t> xyz_;
301 RCP<tMVector_t> vtxWeights_;
302 RCP<tMVector_t> edgWeights_;
304 #ifdef HAVE_EPETRA_DATA_TYPES 305 RCP<const Epetra_Comm> ecomm_;
306 RCP<Epetra_CrsMatrix> eM_;
307 RCP<Epetra_CrsGraph> eG_;
315 void readMatrixMarketFile(
string path,
string testData,
bool distributeInput =
true);
320 void buildCrsMatrix(
int xdim,
int ydim,
int zdim,
string type,
321 bool distributeInput);
327 void readZoltanTestData(
string path,
string testData,
328 bool distributeInput);
331 void getUIChacoGraph(FILE *fptr,
bool haveAssign, FILE *assignFile,
332 string name,
bool distributeInput);
335 void getUIChacoCoords(FILE *fptr,
string name);
341 static const int CHACO_LINE_LENGTH=200;
342 char chaco_line[CHACO_LINE_LENGTH];
347 double chaco_read_val(FILE* infile,
int *end_flag);
348 int chaco_read_int(FILE* infile,
int *end_flag);
349 void chaco_flush_line(FILE*);
352 int chaco_input_graph(FILE *fin,
const char *inname,
int **start,
353 int **adjacency,
int *nvtxs,
int *nVwgts,
354 float **vweights,
int *nEwgts,
float **eweights);
357 int chaco_input_geom(FILE *fingeom,
const char *geomname,
int nvtxs,
358 int *igeom,
double **x,
double **y,
double **z);
361 int chaco_input_assign(FILE *finassign,
const char *assignname,
int nvtxs,
369 void readGeometricGenTestData(
string path,
string testData);
373 ParameterList &geoparams);
378 const string& delimiters =
" \f\n\r\t\v" );
381 const string& delimiters =
" \f\n\r\t\v" );
384 const string& delimiters =
" \f\n\r\t\v" );
388 void readPamgenMeshFile(
string path,
string testData);
389 #ifdef HAVE_ZOLTAN2_PAMGEN 390 void setPamgenAdjacencyGraph();
391 void setPamgenCoordinateMV();
396 const RCP<
const Comm<int> > &c,
397 bool debugInfo,
bool distributeInput):
398 verbose_(debugInfo), tcomm_(c), havePamgenMesh(false),
399 M_(), xM_(), xyz_(), vtxWeights_(), edgWeights_(),
401 ecomm_(), eM_(), eG_(),
403 chaco_offset(0), chaco_break_pnt(CHACO_LINE_LENGTH)
405 bool zoltan1 =
false;
406 string::size_type loc = path.find(
"/zoltan/test/");
407 if (loc != string::npos)
411 readZoltanTestData(path, testData, distributeInput);
413 readMatrixMarketFile(path, testData);
415 #ifdef HAVE_EPETRA_DATA_TYPES 416 ecomm_ = Xpetra::toEpetra(c);
422 const RCP<
const Comm<int> > &c,
424 bool distributeInput):
425 verbose_(debugInfo), tcomm_(c), havePamgenMesh(false),
426 M_(), xM_(), xyz_(), vtxWeights_(), edgWeights_(),
428 ecomm_(), eM_(), eG_(),
430 chaco_offset(0), chaco_break_pnt(CHACO_LINE_LENGTH)
432 if (matrixType.size() == 0){
438 matrixType = string(
"Laplace1D");
440 matrixType = string(
"Laplace2D");
442 matrixType = string(
"Laplace3D");
444 throw std::runtime_error(
"input");
446 if (verbose_ && tcomm_->getRank() == 0)
447 std::cout <<
"UserInputForTests, Matrix type : " << matrixType << std::endl;
450 buildCrsMatrix(x, y, z, matrixType, distributeInput);
452 #ifdef HAVE_EPETRA_DATA_TYPES 453 ecomm_ = Xpetra::toEpetra(c);
458 const RCP<
const Comm<int> > &c):
459 tcomm_(c), havePamgenMesh(false),
460 M_(), xM_(), xyz_(), vtxWeights_(), edgWeights_(),
462 ecomm_(), eM_(), eG_(),
464 chaco_offset(0), chaco_break_pnt(CHACO_LINE_LENGTH)
468 bool distributeInput =
true, debugInfo =
true;
470 if(pList.isParameter(
"distribute input"))
471 distributeInput = pList.get<
bool>(
"distribute input");
473 if(pList.isParameter(
"debug"))
474 debugInfo = pList.get<
bool>(
"debug");
475 this->verbose_ = debugInfo;
477 if(pList.isParameter(
"input file"))
482 if(pList.isParameter(
"input path"))
483 path = pList.get<
string>(
"input path");
485 string testData = pList.get<
string>(
"input file");
491 if(pList.isParameter(
"file type") && pList.get<
string>(
"file type") ==
"Geometric Generator")
493 else if(pList.isParameter(
"file type") && pList.get<
string>(
"file type") ==
"Pamgen")
497 else if(pList.isParameter(
"file type") && pList.get<
string>(
"file type") ==
"Chaco")
501 switch (file_format) {
502 case GEOMGEN: readGeometricGenTestData(path,testData);
break;
503 case PAMGEN: readPamgenMeshFile(path,testData);
break;
504 case CHACO: readZoltanTestData(path, testData, distributeInput);
break;
505 default: readMatrixMarketFile(path, testData, distributeInput);
break;
508 }
else if(pList.isParameter(
"x") || pList.isParameter(
"y") || pList.isParameter(
"z")){
512 if(pList.isParameter(
"x")) x = pList.get<
int>(
"x");
513 if(pList.isParameter(
"y")) y = pList.get<
int>(
"y");
514 if(pList.isParameter(
"z")) z = pList.get<
int>(
"z");
516 string problemType =
"";
517 if(pList.isParameter(
"equation type")) problemType = pList.get<
string>(
"equation type");
519 if (problemType.size() == 0){
525 problemType = string(
"Laplace1D");
527 problemType = string(
"Laplace2D");
529 problemType = string(
"Laplace3D");
531 throw std::runtime_error(
"input");
533 if (verbose_ && tcomm_->getRank() == 0)
534 std::cout <<
"UserInputForTests, Matrix type : " << problemType << std::endl;
538 buildCrsMatrix(x, y, z, problemType, distributeInput);
541 std::cerr <<
"Input file block undefined!" << std::endl;
544 #ifdef HAVE_EPETRA_DATA_TYPES 545 ecomm_ = Xpetra::toEpetra(c);
554 throw std::runtime_error(
"could not read coord file");
571 throw std::runtime_error(
"could not read mtx file");
578 throw std::runtime_error(
"could not read mtx file");
579 return rcp_const_cast<
tcrsGraph_t>(M_->getCrsGraph());
584 RCP<tVector_t> V = rcp(
new tVector_t(M_->getRowMap(), 1));
592 RCP<tMVector_t> mV = rcp(
new tMVector_t(M_->getRowMap(), nvec));
601 throw std::runtime_error(
"could not read mtx file");
608 throw std::runtime_error(
"could not read mtx file");
609 return rcp_const_cast<
xcrsGraph_t>(xM_->getCrsGraph());
623 #ifdef HAVE_EPETRA_DATA_TYPES 624 RCP<Epetra_CrsGraph> UserInputForTests::getUIEpetraCrsGraph()
627 throw std::runtime_error(
"could not read mtx file");
628 RCP<const tcrsGraph_t> tgraph = M_->getCrsGraph();
629 RCP<const Tpetra::Map<zlno_t, zgno_t> > trowMap = tgraph->getRowMap();
630 RCP<const Tpetra::Map<zlno_t, zgno_t> > tcolMap = tgraph->getColMap();
632 int nElts =
static_cast<int>(trowMap->getGlobalNumElements());
633 int nMyElts =
static_cast<int>(trowMap->getNodeNumElements());
634 int base = trowMap->getIndexBase();
635 ArrayView<const int> gids = trowMap->getNodeElementList();
637 Epetra_BlockMap erowMap(nElts, nMyElts,
638 gids.getRawPtr(), 1, base, *ecomm_);
640 Array<int> rowSize(nMyElts);
641 for (
int i=0; i < nMyElts; i++){
642 rowSize[i] =
static_cast<int>(M_->getNumEntriesInLocalRow(i+base));
645 size_t maxRow = M_->getNodeMaxNumRowEntries();
646 Array<int> colGids(maxRow);
647 ArrayView<const int> colLid;
649 eG_ = rcp(
new Epetra_CrsGraph(Copy, erowMap,
650 rowSize.getRawPtr(),
true));
652 for (
int i=0; i < nMyElts; i++){
653 tgraph->getLocalRowView(i+base, colLid);
654 for (
int j=0; j < colLid.size(); j++)
655 colGids[j] = tcolMap->getGlobalElement(colLid[j]);
656 eG_->InsertGlobalIndices(gids[i], rowSize[i], colGids.getRawPtr());
662 RCP<Epetra_CrsMatrix> UserInputForTests::getUIEpetraCrsMatrix()
665 throw std::runtime_error(
"could not read mtx file");
666 RCP<Epetra_CrsGraph> egraph = getUIEpetraCrsGraph();
667 eM_ = rcp(
new Epetra_CrsMatrix(Copy, *egraph));
669 size_t maxRow = M_->getNodeMaxNumRowEntries();
670 int nrows = egraph->NumMyRows();
671 int base = egraph->IndexBase();
672 const Epetra_BlockMap &rowMap = egraph->RowMap();
673 const Epetra_BlockMap &colMap = egraph->ColMap();
674 Array<int> colGid(maxRow);
676 for (
int i=0; i < nrows; i++){
677 ArrayView<const int> colLid;
678 ArrayView<const zscalar_t> nz;
679 M_->getLocalRowView(i+base, colLid, nz);
680 size_t rowSize = colLid.size();
681 int rowGid = rowMap.GID(i+base);
682 for (
size_t j=0; j < rowSize; j++){
683 colGid[j] = colMap.GID(colLid[j]);
685 eM_->InsertGlobalValues(rowGid, (
int)rowSize, nz.getRawPtr(), colGid.getRawPtr());
691 RCP<Epetra_Vector> UserInputForTests::getUIEpetraVector()
693 RCP<Epetra_CrsGraph> egraph = getUIEpetraCrsGraph();
694 RCP<Epetra_Vector> V = rcp(
new Epetra_Vector(egraph->RowMap()));
699 RCP<Epetra_MultiVector> UserInputForTests::getUIEpetraMultiVector(
int nvec)
701 RCP<Epetra_CrsGraph> egraph = getUIEpetraCrsGraph();
702 RCP<Epetra_MultiVector> mV =
703 rcp(
new Epetra_MultiVector(egraph->RowMap(), nvec));
713 this->hasUITpetraCrsMatrix() || \
714 this->hasUITpetraCrsGraph() || \
715 this->hasPamgenMesh();
720 if(input_type ==
"coordinates")
722 else if(input_type ==
"tpetra_vector")
724 else if(input_type ==
"tpetra_multivector")
726 else if(input_type ==
"tpetra_crs_graph")
728 else if(input_type ==
"tpetra_crs_matrix")
730 else if(input_type ==
"xpetra_vector")
732 else if(input_type ==
"xpetra_multivector")
734 else if(input_type ==
"xpetra_crs_graph")
736 else if(input_type ==
"xpetra_crs_matrix")
738 #ifdef HAVE_EPETRA_DATA_TYPES 739 else if(input_type ==
"epetra_vector")
740 return this->hasUIEpetraVector();
741 else if(input_type ==
"epetra_multivector")
742 return this->hasUIEpetraMultiVector();
743 else if(input_type ==
"epetra_crs_graph")
744 return this->hasUIEpetraCrsGraph();
745 else if(input_type ==
"epetra_crs_matrix")
746 return this->hasUIEpetraCrsMatrix();
754 return xyz_.is_null() ?
false :
true;
759 return vtxWeights_.is_null() ?
false :
true;
764 return edgWeights_.is_null() ?
false :
true;
769 return M_.is_null() ?
false :
true;
774 return M_.is_null() ?
false :
true;
789 return M_.is_null() ?
false :
true;
794 return M_.is_null() ?
false :
true;
809 return this->havePamgenMesh;
812 #ifdef HAVE_EPETRA_DATA_TYPES 813 bool UserInputForTests::hasUIEpetraCrsGraph()
815 return M_.is_null() ?
false :
true;
818 bool UserInputForTests::hasUIEpetraCrsMatrix()
820 return hasUIEpetraCrsGraph();
823 bool UserInputForTests::hasUIEpetraVector()
825 return hasUIEpetraCrsGraph();
828 bool UserInputForTests::hasUIEpetraMultiVector()
830 return hasUIEpetraCrsGraph();
836 ArrayView<ArrayRCP<zscalar_t > > data)
841 size_t dim = data.size();
842 for (
size_t i=0; i < dim; i++){
845 throw (std::bad_alloc());
846 data[i] = Teuchos::arcp(tmp, 0, length,
true);
849 zscalar_t scalingFactor = (max-min) / RAND_MAX;
851 for (
size_t i=0; i < dim; i++){
853 for (
zlno_t j=0; j < length; j++)
854 *x++ = min + (
zscalar_t(rand()) * scalingFactor);
860 string UserInputForTests::trim_right_copy(
862 const string& delimiters)
864 return s.substr( 0, s.find_last_not_of( delimiters ) + 1 );
867 string UserInputForTests::trim_left_copy(
869 const string& delimiters)
871 return s.substr( s.find_first_not_of( delimiters ) );
874 string UserInputForTests::trim_copy(
876 const string& delimiters)
878 return trim_left_copy( trim_right_copy( s, delimiters ), delimiters );
881 void UserInputForTests::readGeometricGenTestData(
string path,
885 std::ostringstream fname;
886 fname << path <<
"/" << testData <<
".txt";
889 if (verbose_ && tcomm_->getRank() == 0)
890 std::cout <<
"UserInputForTests, Read: " << fname.str() << std::endl;
892 Teuchos::ParameterList geoparams(
"geo params");
893 readGeoGenParams(fname.str(),geoparams);
898 int coord_dim = gg->getCoordinateDimension();
899 int numWeightsPerCoord = gg->getNumWeights();
900 zlno_t numLocalPoints = gg->getNumLocalCoords();
901 zgno_t numGlobalPoints = gg->getNumGlobalCoords();
905 for(
int i = 0; i < coord_dim; ++i){
906 coords[i] =
new zscalar_t[numLocalPoints];
910 gg->getLocalCoordinatesCopy(coords);
914 if (numWeightsPerCoord) {
916 weight =
new zscalar_t * [numWeightsPerCoord];
917 for(
int i = 0; i < numWeightsPerCoord; ++i){
918 weight[i] =
new zscalar_t[numLocalPoints];
922 gg->getLocalWeightsCopy(weight);
929 RCP<Tpetra::Map<zlno_t, zgno_t, znode_t> > mp =
930 rcp(
new Tpetra::Map<zlno_t, zgno_t, znode_t>(numGlobalPoints, numLocalPoints, 0, this->tcomm_));
933 Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > coordView(coord_dim);
934 for (
int i=0; i < coord_dim; i++){
935 if(numLocalPoints > 0){
936 Teuchos::ArrayView<const zscalar_t> a(coords[i], numLocalPoints);
940 Teuchos::ArrayView<const zscalar_t> a;
946 xyz_ = RCP<tMVector_t>(
new 951 if (numWeightsPerCoord) {
953 Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > weightView(numWeightsPerCoord);
954 for (
int i=0; i < numWeightsPerCoord; i++){
955 if(numLocalPoints > 0){
956 Teuchos::ArrayView<const zscalar_t> a(weight[i], numLocalPoints);
960 Teuchos::ArrayView<const zscalar_t> a;
965 vtxWeights_ = RCP<tMVector_t>(
new tMVector_t(mp, weightView.view(0, numWeightsPerCoord),
966 numWeightsPerCoord));
970 void UserInputForTests::readGeoGenParams(
string paramFileName,
971 ParameterList &geoparams){
975 std::string input =
"";
977 for(
int i = 0; i < 25000; ++i){
982 if(this->tcomm_->getRank() == 0){
984 fstream inParam(paramFileName.c_str());
991 std::string tmp =
"";
992 getline (inParam,tmp);
993 while (!inParam.eof()){
995 tmp = trim_copy(tmp);
1000 getline (inParam,tmp);
1003 for (
size_t i = 0; i < input.size(); ++i){
1011 int size = (int)input.size();
1015 this->tcomm_->broadcast(0,
sizeof(
int), (
char*) &size);
1017 throw "File " + paramFileName +
" cannot be opened.";
1019 this->tcomm_->broadcast(0, size, inp);
1020 istringstream inParam(inp);
1022 getline (inParam,str);
1023 while (!inParam.eof()){
1024 if(str[0] != param_comment){
1025 size_t pos = str.find(
'=');
1026 if(pos == string::npos){
1027 throw "Invalid Line:" + str +
" in parameter file";
1029 string paramname = trim_copy(str.substr(0,pos));
1030 string paramvalue = trim_copy(str.substr(pos + 1));
1031 geoparams.set(paramname, paramvalue);
1033 getline (inParam,str);
1037 void UserInputForTests::readMatrixMarketFile(
string path,
string testData,
bool distributeInput)
1039 std::ostringstream fname;
1040 fname << path <<
"/" << testData <<
".mtx";
1043 if (verbose_ && tcomm_->getRank() == 0)
1044 std::cout <<
"UserInputForTests, Read: " << fname.str() << std::endl;
1051 RCP<tcrsMatrix_t> toMatrix;
1052 RCP<tcrsMatrix_t> fromMatrix;
1055 fromMatrix = Tpetra::MatrixMarket::Reader<tcrsMatrix_t>::readSparseFile(
1056 fname.str(), tcomm_, dnode,
true,
true,
false);
1057 if(!distributeInput)
1059 if (verbose_ && tcomm_->getRank() == 0)
1060 std::cout <<
"Constructing serial distribution of matrix" << std::endl;
1062 RCP<const map_t> fromMap = fromMatrix->getRowMap();
1064 size_t numGlobalCoords = fromMap->getGlobalNumElements();
1065 size_t numLocalCoords = this->tcomm_->getRank() == 0 ? numGlobalCoords : 0;
1066 RCP<const map_t> toMap = rcp(
new map_t(numGlobalCoords,numLocalCoords, 0, tcomm_));
1068 RCP<import_t> importer = rcp(
new import_t(fromMap, toMap));
1070 toMatrix->doImport(*fromMatrix, *importer, Tpetra::INSERT);
1071 toMatrix->fillComplete();
1074 toMatrix = fromMatrix;
1076 }
catch (std::exception &e) {
1077 if (tcomm_->getRank() == 0)
1078 std::cout <<
"UserInputForTests unable to read matrix market file:" 1079 << fname.str() << std::endl;
1083 "UserInputForTests unable to read matrix market file");
1092 fname << path <<
"/" << testData <<
"_coord.mtx";
1094 size_t coordDim = 0, numGlobalCoords = 0;
1095 size_t msg[2]={0,0};
1096 ArrayRCP<ArrayRCP<zscalar_t> > xyz;
1097 std::ifstream coordFile;
1099 if (tcomm_->getRank() == 0){
1102 std::cout <<
"UserInputForTests, Read: " <<
1103 fname.str() << std::endl;
1107 coordFile.open(fname.str().c_str());
1109 catch (std::exception &e){
1120 while (!done && !fail && coordFile.good()){
1121 coordFile.getline(c, 256);
1124 else if (c[0] ==
'%')
1128 std::istringstream s(c);
1129 s >> numGlobalCoords >> coordDim;
1130 if (!s.eof() || numGlobalCoords < 1 || coordDim < 1)
1139 xyz = Teuchos::arcp(
new ArrayRCP<zscalar_t> [coordDim], 0, coordDim);
1141 for (
size_t dim=0; !fail && dim < coordDim; dim++){
1147 xyz[dim] = Teuchos::arcp(tmp, 0, numGlobalCoords);
1149 for (idx=0; !coordFile.eof() && idx < numGlobalCoords; idx++){
1150 coordFile.getline(c, 256);
1151 std::istringstream s(c);
1155 if (idx < numGlobalCoords)
1161 ArrayRCP<zscalar_t> emptyArray;
1162 for (
size_t dim=0; dim < coordDim; dim++)
1163 xyz[dim] = emptyArray;
1176 msg[1] = numGlobalCoords;
1180 Teuchos::broadcast<int, size_t>(*tcomm_, 0, 2, msg);
1183 numGlobalCoords = msg[1];
1189 RCP<const map_t> toMap;
1192 base = M_->getIndexBase();
1193 const RCP<const map_t> &mapM = M_->getRowMap();
1197 if (verbose_ && tcomm_->getRank() == 0)
1199 std::cout <<
"Matrix was null. ";
1200 std::cout <<
"Constructing distribution map for coordinate vector." << std::endl;
1204 if(!distributeInput)
1206 if (verbose_ && tcomm_->getRank() == 0)
1207 std::cout <<
"Constructing serial distribution map for coordinates." << std::endl;
1209 size_t numLocalCoords = this->tcomm_->getRank() == 0 ? numGlobalCoords : 0;
1210 toMap = rcp(
new map_t(numGlobalCoords,numLocalCoords, base, tcomm_));
1212 toMap = rcp(
new map_t(numGlobalCoords, base, tcomm_));
1220 ArrayRCP<ArrayView<const zscalar_t> > coordLists(coordDim);
1222 if (tcomm_->getRank() == 0){
1224 for (
size_t dim=0; dim < coordDim; dim++)
1225 coordLists[dim] = xyz[dim].view(0, numGlobalCoords);
1229 throw std::bad_alloc();
1231 ArrayRCP<const zgno_t> rowIds = Teuchos::arcp(tmp, 0, numGlobalCoords);
1233 zgno_t basePlusNumGlobalCoords = base +
static_cast<zgno_t>(numGlobalCoords);
1234 for (
zgno_t id=base;
id < basePlusNumGlobalCoords;
id++)
1237 RCP<const map_t> fromMap = rcp(
new map_t(numGlobalCoords,
1238 rowIds.view(0, numGlobalCoords), base, tcomm_));
1240 tMVector_t allCoords(fromMap, coordLists.view(0, coordDim), coordDim);
1244 xyz_->doExport(allCoords, exporter, Tpetra::INSERT);
1248 RCP<const map_t> fromMap = rcp(
new map_t(numGlobalCoords,
1249 ArrayView<zgno_t>(), base, tcomm_));
1251 tMVector_t allCoords(fromMap, coordLists.view(0, coordDim), coordDim);
1255 xyz_->doExport(allCoords, exporter, Tpetra::INSERT);
1259 void UserInputForTests::buildCrsMatrix(
int xdim,
int ydim,
int zdim,
1260 string problemType,
bool distributeInput)
1262 #ifdef HAVE_ZOLTAN2_GALERI 1263 Teuchos::CommandLineProcessor tclp;
1264 Galeri::Xpetra::Parameters<zgno_t> params(tclp,
1265 xdim, ydim, zdim, problemType);
1267 RCP<const Tpetra::Map<zlno_t, zgno_t> > map;
1268 if (distributeInput)
1269 map = rcp(
new Tpetra::Map<zlno_t, zgno_t>(params.GetNumGlobalElements(),
1273 size_t nGlobalElements = params.GetNumGlobalElements();
1274 size_t nLocalElements = ((tcomm_->getRank() == 0) ? nGlobalElements : 0);
1275 map = rcp(
new Tpetra::Map<zlno_t, zgno_t>(nGlobalElements, nLocalElements, 0,
1279 if (verbose_ && tcomm_->getRank() == 0){
1281 std::cout <<
"Matrix is " << (distributeInput ?
"" :
"not");
1282 std::cout <<
"distributed." << endl;
1284 std::cout <<
"UserInputForTests, Create matrix with " << problemType;
1285 std::cout <<
" (and " << xdim;
1287 std::cout <<
" x " << ydim <<
" x " << zdim;
1289 std::cout <<
" x" << ydim <<
" x 1";
1291 std::cout <<
"x 1 x 1";
1293 std::cout <<
" mesh)" << std::endl;
1299 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 =
1300 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> >
1301 (params.GetMatrixType(), map, params.GetParameterList());
1302 M_ = Pr->BuildMatrix();
1304 catch (std::exception &e) {
1308 "UserInputForTests Galeri::Xpetra::BuildProblem failed");
1314 if (verbose_ && tcomm_->getRank() == 0)
1316 "UserInputForTests, Implied matrix row coordinates computed" <<
1319 ArrayView<const zgno_t> gids = map->getNodeElementList();
1322 size_t pos = problemType.find(
"2D");
1323 if (pos != string::npos)
1325 else if (problemType ==
string(
"Laplace1D") ||
1326 problemType ==
string(
"Identity"))
1332 for (
int i=0; i < dim; i++){
1335 throw(std::bad_alloc());
1336 coordinates[i] = Teuchos::arcp(c, 0, count,
true);
1340 zscalar_t *x = coordinates[0].getRawPtr();
1341 zscalar_t *y = coordinates[1].getRawPtr();
1342 zscalar_t *z = coordinates[2].getRawPtr();
1343 zgno_t xySize = xdim * ydim;
1344 for (
zlno_t i=0; i < count; i++){
1345 zgno_t iz = gids[i] / xySize;
1346 zgno_t xy = gids[i] - iz*xySize;
1353 zscalar_t *x = coordinates[0].getRawPtr();
1354 zscalar_t *y = coordinates[1].getRawPtr();
1355 for (
zlno_t i=0; i < count; i++){
1361 zscalar_t *x = coordinates[0].getRawPtr();
1362 for (
zlno_t i=0; i < count; i++)
1367 Array<ArrayView<const zscalar_t> > coordView(dim);
1369 for (
int i=0; i < dim; i++)
1370 coordView[i] = coordinates[i].view(0,count);
1372 xyz_ = rcp(
new tMVector_t(map, coordView.view(0, dim), dim));
1374 throw std::runtime_error(
"Galeri input requested but Trilinos is " 1375 "not built with Galeri.");
1379 void UserInputForTests::readZoltanTestData(
string path,
string testData,
1380 bool distributeInput)
1383 int rank = tcomm_->getRank();
1384 FILE *graphFile = NULL;
1385 FILE *coordFile = NULL;
1386 FILE *assignFile = NULL;
1391 std::ostringstream chGraphFileName;
1392 chGraphFileName << path <<
"/" << testData <<
".graph";
1395 std::ostringstream chCoordFileName;
1396 chCoordFileName << path <<
"/" << testData <<
".coords";
1399 std::ostringstream chAssignFileName;
1400 chAssignFileName << path <<
"/" << testData <<
".assign";
1403 graphFile = fopen(chGraphFileName.str().c_str(),
"r");
1407 chGraphFileName.str(
"");
1408 chCoordFileName.str(
"");
1410 chGraphFileName << path <<
"/ch_" << testData <<
"/" << testData <<
".graph";
1411 chCoordFileName << path <<
"/ch_" << testData <<
"/" << testData <<
".coords";
1412 chAssignFileName << path <<
"/ch_" << testData <<
"/" << testData <<
".assign";
1415 graphFile = fopen(chGraphFileName.str().c_str(),
"r");
1418 memset(fileInfo, 0,
sizeof(
int) * 3);
1421 if (verbose_ && tcomm_->getRank() == 0)
1422 std::cout <<
"UserInputForTests, open " <<
1423 chGraphFileName.str () << std::endl;
1425 coordFile = fopen(chCoordFileName.str().c_str(),
"r");
1428 if (verbose_ && tcomm_->getRank() == 0)
1429 std::cout <<
"UserInputForTests, open " <<
1430 chCoordFileName.str () << std::endl;
1433 assignFile = fopen(chAssignFileName.str().c_str(),
"r");
1436 if (verbose_ && tcomm_->getRank() == 0)
1437 std::cout <<
"UserInputForTests, open " <<
1438 chAssignFileName.str () << std::endl;
1441 if (verbose_ && tcomm_->getRank() == 0){
1442 std::cout <<
"UserInputForTests, unable to open file: ";
1443 std::cout << chGraphFileName.str() << std::endl;
1449 Teuchos::broadcast<int, int>(*tcomm_, 0, 3, fileInfo);
1451 bool haveGraph = (fileInfo[0] == 1);
1452 bool haveCoords = (fileInfo[1] == 1);
1453 bool haveAssign = (fileInfo[2] == 1);
1458 getUIChacoGraph(graphFile, haveAssign, assignFile,
1459 testData, distributeInput);
1466 getUIChacoCoords(coordFile, testData);
1475 void UserInputForTests::getUIChacoGraph(FILE *fptr,
bool haveAssign,
1476 FILE *assignFile,
string fname,
1477 bool distributeInput)
1479 int rank = tcomm_->getRank();
1481 int nvtxs=0, nedges=0;
1482 int nVwgts=0, nEwgts=0;
1483 int *start = NULL, *adj = NULL;
1484 float *ewgts = NULL, *vwgts = NULL;
1485 size_t *nzPerRow = NULL;
1486 size_t maxRowLen = 0;
1488 ArrayRCP<const size_t> rowSizes;
1490 bool haveEdges =
true;
1494 memset(graphCounts, 0, 5*
sizeof(
int));
1499 fail = chaco_input_graph(fptr, fname.c_str(), &start, &adj,
1500 &nvtxs, &nVwgts, &vwgts, &nEwgts, &ewgts);
1508 std::cout <<
"UserInputForTests, " << nvtxs <<
" vertices,";
1510 std::cout << start[nvtxs] <<
" edges,";
1512 std::cout <<
"no edges,";
1513 std::cout << nVwgts <<
" vertex weights, ";
1514 std::cout << nEwgts <<
" edge weights" << std::endl;
1521 Teuchos::broadcast<int, int>(*tcomm_, 0, 5, graphCounts);
1522 throw std::runtime_error(
"Unable to read chaco file");
1526 nedges = start[nvtxs];
1528 nzPerRow =
new size_t [nvtxs];
1530 throw std::bad_alloc();
1531 rowSizes = arcp(nzPerRow, 0, nvtxs,
true);
1534 for (
int i=0; i < nvtxs; i++){
1535 nzPerRow[i] = start[i+1] - start[i];
1536 if (nzPerRow[i] > maxRowLen)
1537 maxRowLen = nzPerRow[i];
1541 memset(nzPerRow, 0,
sizeof(
size_t) * nvtxs);
1552 int chbase = adj[0];
1553 for (
int i=1; i < nedges; i++)
1554 if (adj[i] < chbase)
1558 for (
int i=0; i < nedges; i++)
1563 graphCounts[0] = nvtxs;
1564 graphCounts[1] = nedges;
1565 graphCounts[2] = nVwgts;
1566 graphCounts[3] = nEwgts;
1567 graphCounts[4] = (int)maxRowLen;
1570 Teuchos::broadcast<int, int>(*tcomm_, 0, 5, graphCounts);
1572 if (graphCounts[0] == 0)
1573 throw std::runtime_error(
"Unable to read chaco file");
1575 haveEdges = (graphCounts[1] > 0);
1577 RCP<tcrsMatrix_t> fromMatrix;
1578 RCP<const map_t> fromMap;
1583 fromMap = rcp(
new map_t(nvtxs, nvtxs, base, tcomm_));
1586 rcp(
new tcrsMatrix_t(fromMap, rowSizes, Tpetra::StaticProfile));
1591 if (nedges && !edgeIds)
1592 throw std::bad_alloc();
1593 for (
int i=0; i < nedges; i++)
1594 edgeIds[i] = adj[i];
1599 zgno_t *nextId = edgeIds;
1600 Array<zscalar_t> values(maxRowLen, 1.0);
1602 for (
int i=0; i < nvtxs; i++){
1603 if (nzPerRow[i] > 0){
1604 ArrayView<const zgno_t> rowNz(nextId, nzPerRow[i]);
1605 fromMatrix->insertGlobalValues(i, rowNz, values.view(0,nzPerRow[i]));
1606 nextId += nzPerRow[i];
1614 fromMatrix->fillComplete();
1617 nvtxs = graphCounts[0];
1618 nedges = graphCounts[1];
1619 nVwgts = graphCounts[2];
1620 nEwgts = graphCounts[3];
1621 maxRowLen = graphCounts[4];
1625 fromMap = rcp(
new map_t(nvtxs, 0, base, tcomm_));
1628 rcp(
new tcrsMatrix_t(fromMap, rowSizes, Tpetra::StaticProfile));
1630 fromMatrix->fillComplete();
1634 RCP<const map_t> toMap;
1635 RCP<tcrsMatrix_t> toMatrix;
1636 RCP<import_t> importer;
1638 if (distributeInput) {
1641 short *assignments =
new short[nvtxs];
1643 fail = chaco_input_assign(assignFile, fname.c_str(), nvtxs, assignments);
1646 Teuchos::broadcast<int, short>(*tcomm_, 0, nvtxs, assignments);
1649 Teuchos::Array<zgno_t> mine;
1650 for (
int i = 0; i < nvtxs; i++) {
1651 if (assignments[i] == rank)
1656 Teuchos::OrdinalTraits<Tpetra::global_size_t>::invalid();
1657 toMap = rcp(
new map_t(dummy, mine(), base, tcomm_));
1658 delete [] assignments;
1662 toMap = rcp(
new map_t(nvtxs, base, tcomm_));
1667 importer = rcp(
new import_t(fromMap, toMap));
1668 toMatrix->doImport(*fromMatrix, *importer, Tpetra::INSERT);
1669 toMatrix->fillComplete();
1673 toMatrix = fromMatrix;
1680 typedef ArrayRCP<const ArrayView<const zscalar_t> > arrayArray_t;
1684 ArrayRCP<zscalar_t> weightBuf;
1685 ArrayView<const zscalar_t> *wgts =
new ArrayView<const zscalar_t> [nVwgts];
1688 size_t len = nVwgts * nvtxs;
1690 if (!buf)
throw std::bad_alloc();
1691 weightBuf = arcp(buf, 0, len,
true);
1693 for (
int widx=0; widx < nVwgts; widx++){
1694 wgts[widx] = ArrayView<const zscalar_t>(buf, nvtxs);
1695 float *vw = vwgts + widx;
1696 for (
int i=0; i < nvtxs; i++, vw += nVwgts)
1705 arrayArray_t vweights = arcp(wgts, 0, nVwgts,
true);
1707 RCP<tMVector_t> fromVertexWeights =
1708 rcp(
new tMVector_t(fromMap, vweights.view(0, nVwgts), nVwgts));
1710 RCP<tMVector_t> toVertexWeights;
1711 if (distributeInput) {
1712 toVertexWeights = rcp(
new tMVector_t(toMap, nVwgts));
1713 toVertexWeights->doImport(*fromVertexWeights, *importer, Tpetra::INSERT);
1716 toVertexWeights = fromVertexWeights;
1718 vtxWeights_ = toVertexWeights;
1723 if (haveEdges && nEwgts > 0){
1725 ArrayRCP<zscalar_t> weightBuf;
1726 ArrayView<const zscalar_t> *wgts =
new ArrayView<const zscalar_t> [nEwgts];
1728 toMap = rcp(
new map_t(nedges, M_->getNodeNumEntries(), base, tcomm_));
1731 size_t len = nEwgts * nedges;
1733 if (!buf)
throw std::bad_alloc();
1734 weightBuf = arcp(buf, 0, len,
true);
1736 for (
int widx=0; widx < nEwgts; widx++){
1737 wgts[widx] = ArrayView<const zscalar_t>(buf, nedges);
1738 float *ew = ewgts + widx;
1739 for (
int i=0; i < nedges; i++, ew += nEwgts)
1746 fromMap = rcp(
new map_t(nedges, nedges, base, tcomm_));
1749 fromMap = rcp(
new map_t(nedges, 0, base, tcomm_));
1752 arrayArray_t eweights = arcp(wgts, 0, nEwgts,
true);
1754 RCP<tMVector_t> fromEdgeWeights;
1755 RCP<tMVector_t> toEdgeWeights;
1756 RCP<import_t> edgeImporter;
1757 if (distributeInput) {
1759 rcp(
new tMVector_t(fromMap, eweights.view(0, nEwgts), nEwgts));
1760 toEdgeWeights = rcp(
new tMVector_t(toMap, nEwgts));
1761 edgeImporter = rcp(
new import_t(fromMap, toMap));
1762 toEdgeWeights->doImport(*fromEdgeWeights, *edgeImporter, Tpetra::INSERT);
1765 toEdgeWeights = fromEdgeWeights;
1767 edgWeights_ = toEdgeWeights;
1771 void UserInputForTests::getUIChacoCoords(FILE *fptr,
string fname)
1773 int rank = tcomm_->getRank();
1775 double *x=NULL, *y=NULL, *z=NULL;
1778 size_t globalNumVtx = M_->getGlobalNumRows();
1785 fail = chaco_input_geom(fptr, fname.c_str(), (int)globalNumVtx,
1792 std::cout <<
"UserInputForTests, read " << globalNumVtx;
1793 std::cout <<
" " << ndim <<
"-dimensional coordinates." << std::endl;
1797 Teuchos::broadcast<int, int>(*tcomm_, 0, 1, &ndim);
1800 throw std::runtime_error(
"Can't read coordinate file");
1802 ArrayRCP<ArrayRCP<const zscalar_t> > coords(ndim);
1807 for (
int dim=0; dim < ndim; dim++){
1810 throw std::bad_alloc();
1811 coords[dim] = arcp<const zscalar_t>(v, 0, globalNumVtx,
true);
1812 double *val = (dim==0 ? x : (dim==1 ? y : z));
1813 for (
size_t i=0; i < globalNumVtx; i++)
1819 len =
static_cast<zlno_t>(globalNumVtx);
1822 RCP<const map_t> fromMap = rcp(
new map_t(globalNumVtx, len, 0, tcomm_));
1823 RCP<const map_t> toMap = M_->getRowMap();
1824 RCP<import_t> importer = rcp(
new import_t(fromMap, toMap));
1826 Array<ArrayView<const zscalar_t> > coordData;
1827 for (
int dim=0; dim < ndim; dim++)
1828 coordData.push_back(coords[dim].view(0, len));
1830 RCP<tMVector_t> fromCoords =
1831 rcp(
new tMVector_t(fromMap, coordData.view(0, ndim), ndim));
1833 RCP<tMVector_t> toCoords = rcp(
new tMVector_t(toMap, ndim));
1835 toCoords->doImport(*fromCoords, *importer, Tpetra::INSERT);
1844 double UserInputForTests::chaco_read_val(
1860 if (chaco_offset == 0 || chaco_offset >= chaco_break_pnt) {
1861 if (chaco_offset >= chaco_break_pnt) {
1862 length_left = CHACO_LINE_LENGTH - chaco_save_pnt - 1;
1864 ptr = &chaco_line[chaco_save_pnt];
1865 for (i=length_left; i; i--) *ptr2++ = *ptr++;
1866 length = chaco_save_pnt + 1;
1869 length = CHACO_LINE_LENGTH;
1874 ptr2 = fgets(&chaco_line[length_left], length, infile);
1876 if (ptr2 == (
char *) NULL) {
1878 return((
double) 0.0);
1881 if ((chaco_line[CHACO_LINE_LENGTH - 2] !=
'\n') && (chaco_line[CHACO_LINE_LENGTH - 2] !=
'\f')
1882 && (strlen(chaco_line) == CHACO_LINE_LENGTH - 1)){
1884 chaco_break_pnt = CHACO_LINE_LENGTH - 1;
1885 chaco_save_pnt = chaco_break_pnt;
1890 if (chaco_line[chaco_break_pnt] !=
'\0') {
1891 if (isspace((
int)(chaco_line[chaco_break_pnt]))) {
1893 chaco_save_pnt = chaco_break_pnt + 1;
1897 else if (white_seen) {
1904 chaco_break_pnt = CHACO_LINE_LENGTH;
1910 while (isspace((
int)(chaco_line[chaco_offset])) && chaco_offset < CHACO_LINE_LENGTH) chaco_offset++;
1911 if (chaco_line[chaco_offset] ==
'%' || chaco_line[chaco_offset] ==
'#') {
1913 if (chaco_break_pnt < CHACO_LINE_LENGTH) {
1914 chaco_flush_line(infile);
1916 return((
double) 0.0);
1919 ptr = &(chaco_line[chaco_offset]);
1920 val = strtod(ptr, &ptr2);
1925 return((
double) 0.0);
1928 chaco_offset = (int) (ptr2 - chaco_line) /
sizeof(char);
1935 int UserInputForTests::chaco_read_int(
1951 if (chaco_offset == 0 || chaco_offset >= chaco_break_pnt) {
1952 if (chaco_offset >= chaco_break_pnt) {
1953 length_left = CHACO_LINE_LENGTH - chaco_save_pnt - 1;
1955 ptr = &chaco_line[chaco_save_pnt];
1956 for (i=length_left; i; i--) *ptr2++ = *ptr++;
1957 length = chaco_save_pnt + 1;
1960 length = CHACO_LINE_LENGTH;
1965 ptr2 = fgets(&chaco_line[length_left], length, infile);
1967 if (ptr2 == (
char *) NULL) {
1972 if ((chaco_line[CHACO_LINE_LENGTH - 2] !=
'\n') && (chaco_line[CHACO_LINE_LENGTH - 2] !=
'\f')
1973 && (strlen(chaco_line) == CHACO_LINE_LENGTH - 1)){
1975 chaco_break_pnt = CHACO_LINE_LENGTH - 1;
1976 chaco_save_pnt = chaco_break_pnt;
1981 if (chaco_line[chaco_break_pnt] !=
'\0') {
1982 if (isspace((
int)(chaco_line[chaco_break_pnt]))) {
1984 chaco_save_pnt = chaco_break_pnt + 1;
1988 else if (white_seen) {
1995 chaco_break_pnt = CHACO_LINE_LENGTH;
2001 while (isspace((
int)(chaco_line[chaco_offset])) && chaco_offset < CHACO_LINE_LENGTH) chaco_offset++;
2002 if (chaco_line[chaco_offset] ==
'%' || chaco_line[chaco_offset] ==
'#') {
2004 if (chaco_break_pnt < CHACO_LINE_LENGTH) {
2005 chaco_flush_line(infile);
2010 ptr = &(chaco_line[chaco_offset]);
2011 val = (int) strtol(ptr, &ptr2, 10);
2019 chaco_offset = (int) (ptr2 - chaco_line) /
sizeof(char);
2025 void UserInputForTests::chaco_flush_line(
2032 while (c !=
'\n' && c !=
'\f')
2036 int UserInputForTests::chaco_input_graph(
2086 while (end_flag == 1) {
2087 *nvtxs = chaco_read_int(fin, &end_flag);
2091 printf(
"ERROR in graph file `%s':", inname);
2092 printf(
" Invalid number of vertices (%d).\n", *nvtxs);
2097 narcs = chaco_read_int(fin, &end_flag);
2099 printf(
"ERROR in graph file `%s':", inname);
2100 printf(
" Invalid number of expected edges (%d).\n", narcs);
2107 option = chaco_read_int(fin, &end_flag);
2109 using_ewgts = option - 10 * (option / 10);
2111 using_vwgts = option - 10 * (option / 10);
2113 vtxnums = option - 10 * (option / 10);
2116 (*nVwgts) = using_vwgts;
2117 (*nEwgts) = using_ewgts;
2120 if (!end_flag && using_vwgts==1){
2121 j = chaco_read_int(fin, &end_flag);
2122 if (!end_flag) (*nVwgts) = j;
2124 if (!end_flag && using_ewgts==1){
2125 j = chaco_read_int(fin, &end_flag);
2126 if (!end_flag) (*nEwgts) = j;
2131 j = chaco_read_int(fin, &end_flag);
2134 *start = (
int *) malloc((
unsigned) (*nvtxs + 1) *
sizeof(
int));
2136 *adjacency = (
int *) malloc((
unsigned) (2 * narcs + 1) *
sizeof(
int));
2141 *vweights = (
float *) malloc((
unsigned) (*nvtxs) * (*nVwgts) *
sizeof(float));
2146 *eweights = (
float *)
2147 malloc((
unsigned) (2 * narcs + 1) * (*nEwgts) *
sizeof(float));
2151 adjptr = *adjacency;
2160 while ((using_vwgts || vtxnums || narcs) && end_flag != -1) {
2165 j = chaco_read_int(fin, &end_flag);
2167 if (vertex == *nvtxs)
2169 printf(
"ERROR in graph file `%s':", inname);
2170 printf(
" no vertex number in line %d.\n", line_num);
2174 if (j != vertex && j != vertex + 1) {
2175 printf(
"ERROR in graph file `%s':", inname);
2176 printf(
" out-of-order vertex number in line %d.\n", line_num);
2190 if (vertex > *nvtxs)
2194 if (using_vwgts && new_vertex) {
2195 for (j=0; j<(*nVwgts); j++){
2196 weight = chaco_read_val(fin, &end_flag);
2198 printf(
"ERROR in graph file `%s':", inname);
2199 printf(
" not enough weights for vertex %d.\n", vertex);
2203 (*vweights)[(vertex-1)*(*nVwgts)+j] = weight;
2210 neighbor = chaco_read_int(fin, &end_flag);
2216 for (j=0; j<(*nEwgts); j++){
2217 eweight = chaco_read_val(fin, &end_flag);
2220 printf(
"ERROR in graph file `%s':", inname);
2221 printf(
" not enough weights for edge (%d,%d).\n", vertex, neighbor);
2234 if (++nedges > 2*narcs) {
2235 printf(
"ERROR in graph file `%s':", inname);
2236 printf(
" at least %d adjacencies entered, but nedges = %d\n",
2241 *adjptr++ = neighbor;
2246 neighbor = chaco_read_int(fin, &end_flag);
2250 (*start)[vertex] = sum_edges;
2255 while (!flag && end_flag != -1) {
2256 chaco_read_int(fin, &end_flag);
2261 (*start)[*nvtxs] = sum_edges;
2269 if (*adjacency != NULL)
2272 if (*eweights != NULL)
2281 if (*adjacency != NULL)
2283 if (*vweights != NULL)
2285 if (*eweights != NULL)
2293 return (error_flag);
2297 int UserInputForTests::chaco_input_geom(
2299 const char *geomname,
2307 double xc, yc, zc =0;
2315 *x = *y = *z = NULL;
2318 while (end_flag == 1) {
2319 xc = chaco_read_val(fingeom, &end_flag);
2323 if (end_flag == -1) {
2324 printf(
"No values found in geometry file `%s'\n", geomname);
2330 yc = chaco_read_val(fingeom, &end_flag);
2331 if (end_flag == 0) {
2333 zc = chaco_read_val(fingeom, &end_flag);
2334 if (end_flag == 0) {
2336 chaco_read_val(fingeom, &end_flag);
2338 printf(
"Too many values on input line of geometry file `%s'\n",
2341 printf(
" Maximum dimensionality is 3\n");
2350 *x = (
double *) malloc((
unsigned) nvtxs *
sizeof(double));
2353 *y = (
double *) malloc((
unsigned) nvtxs *
sizeof(double));
2357 *z = (
double *) malloc((
unsigned) nvtxs *
sizeof(double));
2361 for (nread = 1; nread < nvtxs; nread++) {
2364 i = fscanf(fingeom,
"%lf", &((*x)[nread]));
2366 else if (ndims == 2) {
2367 i = fscanf(fingeom,
"%lf%lf", &((*x)[nread]), &((*y)[nread]));
2369 else if (ndims == 3) {
2370 i = fscanf(fingeom,
"%lf%lf%lf", &((*x)[nread]), &((*y)[nread]),
2375 printf(
"Too few lines of values in geometry file; nvtxs=%d, but only %d read\n",
2380 else if (i != ndims) {
2381 printf(
"Wrong number of values in line %d of geometry file `%s'\n",
2382 line_num, geomname);
2391 while (!flag && end_flag != -1) {
2392 chaco_read_val(fingeom, &end_flag);
2404 int UserInputForTests::chaco_input_assign(
2406 const char *inassignname,
2417 while (end_flag == 1) {
2418 assignment[0] = chaco_read_int(finassign, &end_flag);
2421 if (assignment[0] < 0) {
2422 printf(
"ERROR: Entry %d in assignment file `%s' less than zero (%d)\n",
2423 1, inassignname, assignment[0]);
2428 if (end_flag == -1) {
2429 printf(
"ERROR: No values found in assignment file `%s'\n", inassignname);
2435 if (assignment[0] > nvtxs)
2436 flag = assignment[1];
2437 for (i = 1; i < nvtxs; i++) {
2438 j = fscanf(finassign,
"%hd", &(assignment[i]));
2440 printf(
"ERROR: Too few values in assignment file `%s'.\n", inassignname);
2444 if (assignment[i] < 0) {
2445 printf(
"ERROR: Entry %d in assignment file `%s' less than zero (%d)\n",
2446 i+1, inassignname, assignment[i]);
2450 if (assignment[i] > nvtxs) {
2451 if (assignment[i] > flag)
2452 flag = assignment[i];
2457 printf(
"WARNING: Possible error in assignment file `%s'\n", inassignname);
2458 printf(
" More assignment sets (%d) than vertices (%d)\n", flag, nvtxs);
2464 while (!flag && end_flag != -1) {
2465 chaco_read_int(finassign, &end_flag);
2470 printf(
"WARNING: Possible error in assignment file `%s'\n", inassignname);
2471 printf(
" Numerical data found after expected end of file\n");
2479 void UserInputForTests::readPamgenMeshFile(
string path,
string testData)
2481 #ifdef HAVE_ZOLTAN2_PAMGEN 2482 int rank = this->tcomm_->getRank();
2483 if (verbose_ && tcomm_->getRank() == 0)
2484 std::cout <<
"UserInputForTestsBD::readPamgenFile, Read: " << testData << std::endl;
2491 std::ostringstream meshFileName;
2492 meshFileName << path <<
"/" << testData <<
".pmgen";
2495 file.open(meshFileName.str(), ios::in);
2499 if(verbose_ && tcomm_->getRank() == 0)
2501 std::cout <<
"Unable to open pamgen mesh: ";
2502 std::cout << meshFileName.str();
2503 std::cout <<
"\nPlease check file path and name." << std::endl;
2509 file.seekg (0,file.end);
2516 while(std::getline(file,line))
2518 if( line.find(
"nz") != std::string::npos ||
2519 line.find(
"nphi") != std::string::npos)
2527 file.seekg(0, ios::beg);
2532 this->tcomm_->broadcast(0,
sizeof(
int), (
char *)&dimension);
2533 this->tcomm_->broadcast(0,
sizeof(
size_t),(
char *)&len);
2534 this->tcomm_->barrier();
2537 if(verbose_ && tcomm_->getRank() == 0)
2538 std::cout <<
"Pamgen Mesh file size == 0, exiting UserInputForTests early." << endl;
2542 char * file_data =
new char[len];
2543 file_data[len] =
'\0';
2545 file.read(file_data,len);
2549 this->tcomm_->broadcast(0,(
int)len,file_data);
2550 this->tcomm_->barrier();
2554 this->pamgen_mesh = rcp(
new PamgenMesh);
2555 this->havePamgenMesh =
true;
2556 pamgen_mesh->createMesh(file_data,dimension,this->tcomm_);
2559 pamgen_mesh->storeMesh();
2560 this->tcomm_->barrier();
2563 this->setPamgenCoordinateMV();
2566 this->setPamgenAdjacencyGraph();
2568 this->tcomm_->barrier();
2569 if(rank == 0) file.close();
2570 delete [] file_data;
2572 throw std::runtime_error(
"Pamgen requested but Trilinos " 2573 "not built with Pamgen");
2577 #ifdef HAVE_ZOLTAN2_PAMGEN 2578 void UserInputForTests::setPamgenCoordinateMV()
2580 int dimension = pamgen_mesh->num_dim;
2584 zgno_t numelements = pamgen_mesh->num_elem;
2585 zgno_t numGlobalElements = pamgen_mesh->num_elems_global;
2588 for(
int i = 0; i < dimension; ++i){
2589 elem_coords[i] =
new zscalar_t[numelements];
2590 memcpy(elem_coords[i],&pamgen_mesh->element_coord[i*numelements],
sizeof(
double) * numelements);
2594 typedef Tpetra::Map<zlno_t, zgno_t, znode_t>
map_t;
2595 RCP<Tpetra::Map<zlno_t, zgno_t, znode_t> > mp;
2599 Array<zgno_t> elementList(numelements);
2600 for (Array<zgno_t>::size_type k = 0; k < numelements; ++k) {
2601 elementList[k] = pamgen_mesh->element_order_map[k];
2604 mp = rcp (
new map_t (numGlobalElements, elementList, 0, this->tcomm_));
2608 Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > coordView(dimension);
2609 for (
int i = 0; i < dimension; i++){
2610 if(numelements > 0){
2611 Teuchos::ArrayView<const zscalar_t> a(elem_coords[i], numelements);
2615 Teuchos::ArrayView<const zscalar_t> a;
2621 xyz_ = RCP<tMVector_t>(
new 2626 void UserInputForTests::setPamgenAdjacencyGraph()
2634 typedef Tpetra::Map<zlno_t, zgno_t, znode_t>
map_t;
2637 size_t local_nodes = (size_t)this->pamgen_mesh->num_nodes;
2638 size_t local_els = (
size_t)this->pamgen_mesh->num_elem;
2640 size_t global_els = (size_t)this->pamgen_mesh->num_elems_global;
2641 size_t global_nodes = (
size_t)this->pamgen_mesh->num_nodes_global;
2645 RCP<const map_t> rowMap = rcp(
new map_t(global_els,0,this->tcomm_));
2646 RCP<const map_t> rangeMap = rowMap;
2649 RCP<const map_t> domainMap = rcp(
new map_t(global_nodes,0,this->tcomm_));
2652 Teuchos::RCP<tcrsMatrix_t> C = rcp(
new tcrsMatrix_t(rowMap,0));
2655 Array<zgno_t> g_el_ids(local_els);
2656 for (
size_t k = 0; k < local_els; ++k) {
2657 g_el_ids[k] = pamgen_mesh->global_element_numbers[k]-1;
2660 Array<zgno_t> g_node_ids(local_nodes);
2661 for (
size_t k = 0; k < local_nodes; ++k) {
2662 g_node_ids[k] = pamgen_mesh->global_node_numbers[k]-1;
2665 int blks = this->pamgen_mesh->num_elem_blk;
2671 for(
int i = 0; i < blks; i++)
2673 int el_per_block = this->pamgen_mesh->elements[i];
2674 int nodes_per_el = this->pamgen_mesh->nodes_per_element[i];
2675 int * connect = this->pamgen_mesh->elmt_node_linkage[i];
2677 for(
int j = 0; j < el_per_block; j++)
2679 const zgno_t gid =
static_cast<zgno_t>(g_el_ids[el_no]);
2680 for(
int k = 0; k < nodes_per_el; k++)
2682 int g_node_i = g_node_ids[connect[j*nodes_per_el+k]-1];
2683 C->insertGlobalValues(gid,
2684 Teuchos::tuple<zgno_t>(g_node_i),
2685 Teuchos::tuple<zscalar_t>(one));
2691 C->fillComplete(domainMap, rangeMap);
2697 Tpetra::MatrixMatrix::Multiply(*C,
false, *C,
true, *A);
2705 for(
zgno_t gid : rowMap->getNodeElementList())
2707 size_t numEntriesInRow = A->getNumEntriesInGlobalRow (gid);
2708 Array<zscalar_t> rowvals (numEntriesInRow);
2709 Array<zgno_t> rowinds (numEntriesInRow);
2712 Array<zscalar_t> mod_rowvals;
2713 Array<zgno_t> mod_rowinds;
2714 A->getGlobalRowCopy (gid, rowinds (), rowvals (), numEntriesInRow);
2715 for (
size_t i = 0; i < numEntriesInRow; i++) {
2718 if (rowvals[i] >= 1)
2720 mod_rowvals.push_back(one);
2721 mod_rowinds.push_back(rowinds[i]);
2724 this->M_->insertGlobalValues(gid, mod_rowinds, mod_rowvals);
2727 this->M_->fillComplete();
#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.
common code used by tests
#define HAVE_EPETRA_DATA_TYPES
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")
static const std::string fail
string trim_left_copy(const string &s, const string &delimiters=" \f\n\r\t\v")