9#include <fei_macros.hpp>
11#include <test_utils/snl_fei_tester.hpp>
13#include <fei_LinearSystemCore.hpp>
14#include <fei_ArrayUtils.hpp>
15#include <test_utils/LibraryFactory.hpp>
17#include <fei_base.hpp>
20#include <FETI_DP_FiniteElementData.h>
23#include <test_utils/DataReader.hpp>
24#include <test_utils/SolnCheck.hpp>
27#define fei_file "snl_fei_tester.cpp"
28#include <fei_ErrMacros.hpp>
32 MPI_Comm comm,
int localProc,
int numProcs)
52snl_fei_tester::~snl_fei_tester()
59int snl_fei_tester::testInitialization()
61 if (factory_.
get() == NULL) {
65 catch (std::runtime_error& exc) {
69 if (factory_.
get() == NULL) {
74 std::vector<std::string> stdstrings;
90 defineFieldsAndIDTypes();
97 CHK_ERR( initElemBlocks() );
99 CHK_ERR( initConstraints() );
102 for(i=0; i<data_->numSharedNodeSets_; ++i) {
103 CommNodeSet& nodeSet = data_->sharedNodeSets_[i];
106 idTypes_[nodeTypeOffset_],
108 nodeSet.procsPerNode_,
118void snl_fei_tester::dumpMatrixFiles()
120 FEI_OSTRINGSTREAM osstr;
121 osstr <<
"A_" << A_->
typeName() <<
".np"<<numProcs_;
122 std::string str = osstr.str();
127void snl_fei_tester::setParameter(
const char* param)
129 std::vector<std::string> stdstrings;
142int snl_fei_tester::testLoading()
150 matrixGraph_->
setIndicesMode(fei::MatrixGraph::POINT_ENTRY_GRAPH);
152 CHK_ERR( linSys_->
parameters(data_->numParams_, data_->paramStrings_) );
154 std::vector<std::string> stdstrings;
169 CHK_ERR( loadElemBlocks() );
171 CHK_ERR( loadConstraints() );
174 for(i=0; i<data_->numBCNodeSets_; ++i) {
175 BCNodeSet& bcSet = data_->bcNodeSets_[i];
176 int fieldSize = data_->getFieldSize(bcSet.fieldID_);
183 idTypes_[nodeTypeOffset_],
185 bcSet.offsetsIntoField_,
186 bcSet.prescribed_values_) );
195int snl_fei_tester::testSolve()
199 std::vector<std::string> stdstrings;
204 int status, itersTaken = 0;
205 CHK_ERR( solver->solve(linSys_.
get(),
207 paramset, itersTaken, status) );
215int snl_fei_tester::testCheckResult()
217 CHK_ERR( save_block_node_soln(*data_, x_.
get(), data_->solnFileName_.c_str(),
218 numProcs_, localProc_, 1));
220 CHK_ERR( save_block_elem_soln(*data_, x_.
get(), data_->solnFileName_.c_str(),
221 numProcs_, localProc_, 1));
223 CHK_ERR( save_multiplier_soln(*data_, x_.
get(), data_->solnFileName_.c_str(),
224 numProcs_, localProc_, 1));
226 int err = SolnCheck::checkSolution(localProc_, numProcs_, data_->solnFileName_.c_str(),
227 data_->checkFileName_.c_str(),
"node", 1);
229 err += SolnCheck::checkSolution(localProc_, numProcs_, data_->solnFileName_.c_str(),
230 data_->checkFileName_.c_str(),
"elem", 1);
232 err += SolnCheck::checkSolution(localProc_, numProcs_, data_->solnFileName_.c_str(),
233 data_->checkFileName_.c_str(),
"mult", 1);
236 if (MPI_SUCCESS != MPI_Allreduce(&err, &globalErr, 1, MPI_INT, MPI_SUM,
239 if (globalErr != 0)
return(-1);
244void snl_fei_tester::defineFieldsAndIDTypes()
246 vecSpace_->
defineFields(data_->numFields_, data_->fieldIDs_, data_->fieldSizes_);
249 idTypes_.push_back(0);
252 idTypes_.push_back(1);
255 idTypes_.push_back(2);
260 constraintTypeOffset_ = 1;
265int snl_fei_tester::initElemBlocks()
267 for(
int i=0; i<data_->numElemBlocks_; ++i) {
268 ElemBlock& eb = data_->elemBlocks_[i];
271 definePattern(eb, patternID);
277 for(
int j=0; j<eb.numElements_; ++j) {
278 std::vector<int> conn(eb.numNodesPerElement_);
279 for(
int ii=0; ii<eb.numNodesPerElement_; ++ii) {
280 conn[ii] = eb.elemConn_[j][ii];
293int snl_fei_tester::loadElemBlocks()
296 for(i=0; i<data_->numElemBlocks_; ++i) {
297 ElemBlock& eb = data_->elemBlocks_[i];
299 if (eb.numElements_ < 1) {
305 std::vector<int> indices(numIndices);
307 for(
int j=0; j<eb.numElements_; ++j) {
314 if (numIndices != checkNum) {
318 CHK_ERR( A_->
sumIn(eb.blockID_, eb.elemIDs_[j],
321 CHK_ERR( b_->
sumIn(numIndices, &indices[0],
322 eb.elemLoad_[j], 0) );
330int snl_fei_tester::initConstraints()
332 std::vector<int> idTypes;
333 int constraintID = localProc_*100000;
335 for(i=0; i<data_->numCRMultSets_; ++i) {
336 CRSet& crSet = data_->crMultSets_[i];
338 for(
int j=0; j<1; ++j) {
339 idTypes.assign(crSet.
numNodes_, idTypes_[nodeTypeOffset_]);
341 crSet.
crID_ = constraintID++;
342 int constraintIDType = idTypes_[constraintTypeOffset_];
352 for(i=0; i<data_->numCRPenSets_; ++i) {
353 CRSet& crSet = data_->crPenSets_[i];
355 for(
int j=0; j<1; ++j) {
356 idTypes.assign(crSet.
numNodes_, idTypes_[nodeTypeOffset_]);
358 crSet.
crID_ = constraintID++;
359 int constraintIDType = idTypes_[constraintTypeOffset_];
369 std::map<int,int> fieldDB;
370 for(i=0; i<data_->numFields_; ++i) {
371 fieldDB.insert(std::pair<int,int>(data_->fieldIDs_[i], data_->fieldSizes_[i]));
374 std::vector<int> nodeIDs;
375 std::vector<int> fieldIDs;
376 std::vector<double> weights;
378 for(i=0; i<data_->numSlaveVars_; i++) {
380 CRSet& crSet = data_->slaveVars_[i];
388 nodeIDs[ii+1] = crSet.nodeIDs_[0][ii];
389 fieldIDs.push_back(crSet.fieldIDs_[ii]);
392 idTypes.assign(crSet.
numNodes_+1, idTypes_[nodeTypeOffset_]);
396 for(ii=0; ii<fieldSize; ++ii) weights.push_back(0.0);
400 fieldSize = fieldDB[crSet.fieldIDs_[ii]];
401 for(
int jj=0; jj<fieldSize; ++jj) {
402 weights.push_back(crSet.weights_[offset++]);
420int snl_fei_tester::loadConstraints()
423 for(i=0; i<data_->numCRMultSets_; ++i) {
424 CRSet& crSet = data_->crMultSets_[i];
426 for(
int j=0; j<1; ++j) {
433 for(i=0; i<data_->numCRPenSets_; ++i) {
434 CRSet& crSet = data_->crPenSets_[i];
436 for(
int j=0; j<1; ++j) {
448void snl_fei_tester::definePattern(ElemBlock& eb,
int& patternID)
450 int i, j, numIDTypes = 1;
451 numIDTypes += eb.numElemDOF_>0 ? 1 : 0;
454 std::vector<int> nodalFieldIDs;
455 std::vector<int> flatFieldIDsArray;
456 for(i=0; i<eb.numNodesPerElement_; ++i) {
457 for(j=0; j<eb.numFieldsPerNode_[i]; ++j) {
459 flatFieldIDsArray.push_back(eb.nodalFieldIDs_[i][j]);
463 patternID = numPatterns_++;
465 if (numIDTypes == 1 && nodalFieldIDs.size() == 1) {
467 patternID = matrixGraph_->
definePattern(eb.numNodesPerElement_,
468 idTypes_[nodeTypeOffset_],
471 else if (numIDTypes == 1) {
472 std::vector<int> numFieldsPerID(eb.numNodesPerElement_);
474 patternID = matrixGraph_->
definePattern(eb.numNodesPerElement_,
475 idTypes_[nodeTypeOffset_],
476 eb.numFieldsPerNode_,
477 &flatFieldIDsArray[0]);
480 std::vector<int> idTypes(eb.numNodesPerElement_+1, idTypes_[nodeTypeOffset_]);
481 idTypes[idTypes.size()-1] = idTypes_[elemTypeOffset_];
482 std::vector<int> numFieldsPerID(idTypes.size());
483 std::vector<int> fieldIDs;
484 for(i=0; i<eb.numNodesPerElement_; ++i) {
485 numFieldsPerID[i] = eb.numFieldsPerNode_[i];
486 for(j=0; j<eb.numFieldsPerNode_[i]; ++j) {
487 fieldIDs.push_back(eb.nodalFieldIDs_[i][j]);
490 numFieldsPerID[idTypes.size()-1] = eb.numElemDOF_;
491 for(i=0; i<eb.numElemDOF_; ++i) {
492 fieldIDs.push_back(eb.elemDOFFieldIDs_[i]);
503int snl_fei_tester::save_block_node_soln(DataReader& data,
fei::Vector* vec,
504 const char* solnFileName,
int numProcs,
505 int localProc,
int solveCounter)
511 int* nodeList =
new int[numLocalNodes];
515 numLocalNodes, nodeList, checkNum);
520 FEI_OSTRINGSTREAM fileName;
522 std::string str = fileName.str();
523 FEI_OFSTREAM outfile(str.c_str());
525 if (!outfile || outfile.bad()) {
526 fei::console_out() <<
"ERROR opening solution output file " << fileName.str() << FEI_ENDL;
530 outfile.setf(IOS_SCIENTIFIC, IOS_FLOATFIELD);
532 std::vector<double> solnData;
533 std::vector<int> fieldList;
537 for(
int i=0; i<numLocalNodes; i++) {
538 int idType = idTypes_[nodeTypeOffset_];
539 int ID = nodeList[i];
543 solnData.resize(numDOF);
544 vecSpace_->
getFields(idType, ID, fieldList);
546 outfile << ID <<
" " << numDOF << FEI_ENDL;
547 for(
int j=0; j<numFields; ++j) {
549 totalSize += fieldSize;
552 1, &ID, &solnData[0]) );
554 for(
int k=0; k<fieldSize; ++k) {
555 outfile << solnData[k] <<
" ";
561 FEI_COUT <<
"save-node-soln: wrote " << totalSize <<
" entries for " << numLocalNodes <<
" nodes to " << str << FEI_ENDL;
570int snl_fei_tester::save_block_elem_soln(DataReader& data,
fei::Vector* vec,
571 const char* solnFileName,
573 int localProc,
int solveCounter)
579 int* elemList =
new int[numLocalElems];
583 numLocalElems, elemList, checkNum);
588 FEI_OSTRINGSTREAM fileName;
590 std::string str = fileName.str();
591 FEI_OFSTREAM outfile(str.c_str());
593 if (!outfile || outfile.bad()) {
594 fei::console_out() <<
"ERROR opening solution output file " << fileName.str() << FEI_ENDL;
598 std::vector<double> solnData;
599 std::vector<int> fieldList;
601 for(
int i=0; i<numLocalElems; i++) {
602 int idType = idTypes_[elemTypeOffset_];
603 int ID = elemList[i];
607 solnData.resize(numDOF);
608 vecSpace_->
getFields(idType, ID, fieldList);
610 outfile << ID <<
" " << numDOF << FEI_ENDL;
611 for(
int j=0; j<numFields; ++j) {
615 1, &ID, &solnData[0]) );
617 for(
int k=0; k<fieldSize; ++k) {
618 outfile << solnData[k] <<
" ";
631int snl_fei_tester::save_multiplier_soln(DataReader& data,
fei::Vector* vec,
632 const char* solnFileName,
633 int numProcs,
int localProc,
640 int* globalNumCRs =
new int[
numProcs];
642 if (MPI_Allgather(&numLocalCRs, 1, MPI_INT, globalNumCRs, 1, MPI_INT,
643 comm_) != MPI_SUCCESS) {
648 int localCRStart = 0;
650 for(
int p=0; p<
localProc; p++) localCRStart += globalNumCRs[p];
653 delete [] globalNumCRs;
655 std::vector<int> crList(numLocalCRs);
659 idTypes_[constraintTypeOffset_], numLocalCRs,
660 numLocalCRs ? &crList[0] : 0, checkNum);
665 FEI_OSTRINGSTREAM fileName;
667 std::string str = fileName.str();
668 FEI_OFSTREAM outfile(str.c_str());
670 if (!outfile || outfile.bad()) {
671 fei::console_out() <<
"ERROR opening solution output file " << fileName.str() << FEI_ENDL;
675 std::vector<double> solnData;
676 std::vector<int> fieldList;
678 for(
int i=0; i<numLocalCRs; i++) {
679 int idType = idTypes_[constraintTypeOffset_];
684 outfile << localCRStart++ <<
" " << 1 << FEI_ENDL;
685 for(
int j=0; j<1; ++j) {
686 int globalIndex = -1;
689 CHK_ERR( vec->
copyOut(1, &globalIndex, &solnData[0]) );
691 for(
int k=0; k<1; ++k) {
692 outfile << solnData[k] <<
" ";
virtual void parameters(const fei::ParameterSet ¶mset)
virtual fei::SharedPtr< fei::LinearSystem > createLinearSystem(fei::SharedPtr< fei::MatrixGraph > &matrixGraph)
virtual void setRHS(fei::SharedPtr< fei::Vector > &rhs)
virtual int loadEssentialBCs(int numIDs, const int *IDs, int idType, int fieldID, int offsetIntoField, const double *prescribedValues)
virtual void setSolutionVector(fei::SharedPtr< fei::Vector > &soln)
virtual int loadPenaltyConstraint(int constraintID, const double *weights, double penaltyValue, double rhsValue)=0
virtual int loadLagrangeConstraint(int constraintID, const double *weights, double rhsValue)=0
virtual int parameters(int numParams, const char *const *paramStrings)=0
virtual void setMatrix(fei::SharedPtr< fei::Matrix > &matrix)
virtual int loadComplete(bool applyBCs=true, bool globalAssemble=true)=0
virtual fei::SharedPtr< fei::MatrixGraph > createMatrixGraph(fei::SharedPtr< fei::VectorSpace > rowSpace, fei::SharedPtr< fei::VectorSpace > columnSpace, const char *name)=0
virtual int definePattern(int numIDs, int idType)=0
virtual int initPenaltyConstraint(int constraintID, int constraintIDType, int numIDs, const int *idTypes, const int *IDs, const int *fieldIDs)=0
virtual int initComplete()=0
virtual void setIndicesMode(int mode)=0
virtual int createSlaveMatrices()=0
virtual int initConnectivity(int blockID, int connectivityID, const int *connectedIdentifiers)=0
virtual int getConnectivityNumIndices(int blockID) const =0
virtual int initSlaveConstraint(int numIDs, const int *idTypes, const int *IDs, const int *fieldIDs, int offsetOfSlave, int offsetIntoSlaveField, const double *weights, double rhsValue)=0
virtual int initLagrangeConstraint(int constraintID, int constraintIDType, int numIDs, const int *idTypes, const int *IDs, const int *fieldIDs)=0
virtual int getConnectivityIndices(int blockID, int connectivityID, int indicesAllocLen, int *indices, int &numIndices)=0
virtual void setParameters(const fei::ParameterSet ¶ms)=0
virtual int initConnectivityBlock(int blockID, int numConnectivityLists, int patternID, bool diagonal=false)=0
virtual fei::SharedPtr< fei::Matrix > createMatrix(fei::SharedPtr< fei::MatrixGraph > matrixGraph)=0
virtual int writeToFile(const char *filename, bool matrixMarketFormat=true)=0
virtual const char * typeName()=0
virtual int parameters(const fei::ParameterSet ¶mset)=0
virtual int sumIn(int numRows, const int *rows, int numCols, const int *cols, const double *const *values, int format=0)=0
virtual int putScalar(double scalar)=0
void add(const Param ¶m, bool maintain_unique_keys=true)
virtual fei::SharedPtr< fei::Solver > createSolver(const char *name=0)=0
virtual fei::SharedPtr< VectorSpace > createVectorSpace(MPI_Comm, const char *name)
unsigned getFieldSize(int fieldID)
int getNumOwnedAndSharedIDs(int idType)
int initSharedIDs(int numShared, int idType, const int *sharedIDs, const int *numSharingProcsPerID, const int *sharingProcs)
void setParameters(const fei::ParameterSet ¶mset)
void getFields(std::vector< int > &fieldIDs)
int getOwnedAndSharedIDs(int idtype, int lenList, int *IDs, int &numOwnedAndSharedIDs)
int getNumDegreesOfFreedom(int idType, int ID)
int getGlobalIndex(int idType, int ID, int fieldID, int fieldOffset, int whichComponentOfField, int &globalIndex)
void defineFields(int numFields, const int *fieldIDs, const int *fieldSizes, const int *fieldTypes=NULL)
void defineIDTypes(int numIDTypes, const int *idTypes)
virtual fei::SharedPtr< fei::Vector > createVector(fei::SharedPtr< fei::VectorSpace > vecSpace, int numVectors=1)=0
virtual int copyOut(int numValues, const int *indices, double *values, int vectorIndex=0) const =0
virtual int putScalar(double scalar)=0
virtual int scatterToOverlap()=0
virtual int sumIn(int numValues, const int *indices, const double *values, int vectorIndex=0)=0
virtual int copyOutFieldData(int fieldID, int idType, int numIDs, const int *IDs, double *data, int vectorIndex=0)=0
void parse_strings(std::vector< std::string > &stdstrings, const char *separator_string, fei::ParameterSet ¶mset)
void char_ptrs_to_strings(int numStrings, const char *const *charstrings, std::vector< std::string > &stdstrings)
int localProc(MPI_Comm comm)
fei::SharedPtr< fei::Factory > create_fei_Factory(MPI_Comm comm, const char *libraryName)
std::ostream & console_out()
int sortedListInsert(const T &item, std::vector< T > &list)
int numProcs(MPI_Comm comm)