FEI Version of the Day
snl_fei_tester.cpp
00001 /*
00002 // @HEADER
00003 // ************************************************************************
00004 //             FEI: Finite Element Interface to Linear Solvers
00005 //                  Copyright (2005) Sandia Corporation.
00006 //
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, the
00008 // U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Alan Williams (william@sandia.gov) 
00038 //
00039 // ************************************************************************
00040 // @HEADER
00041 */
00042 
00043 #include <fei_macros.hpp>
00044 
00045 #include <test_utils/snl_fei_tester.hpp>
00046 
00047 #include <fei_LinearSystemCore.hpp>
00048 #include <fei_ArrayUtils.hpp>
00049 #include <test_utils/LibraryFactory.hpp>
00050 
00051 #include <fei_base.hpp>
00052 
00053 #ifdef HAVE_FEI_FETI
00054 #include <FETI_DP_FiniteElementData.h>
00055 #endif
00056 
00057 #include <test_utils/DataReader.hpp>
00058 #include <test_utils/SolnCheck.hpp>
00059 
00060 #undef fei_file
00061 #define fei_file "snl_fei_tester.cpp"
00062 #include <fei_ErrMacros.hpp>
00063 
00064 //----------------------------------------------------------------------------
00065 snl_fei_tester::snl_fei_tester(fei::SharedPtr<DataReader> data_reader,
00066                                MPI_Comm comm, int localProc, int numProcs)
00067   : comm_(comm),
00068     factory_(),
00069     vecSpace_(),
00070     matrixGraph_(),
00071     A_(),
00072     x_(),
00073     b_(),
00074     linSys_(NULL),
00075     linSysCore_(NULL),
00076     feData_(NULL),
00077     data_(data_reader),
00078     idTypes_(),
00079     numPatterns_(0),
00080     localProc_(localProc),
00081     numProcs_(numProcs)
00082 {
00083 }
00084 
00085 //----------------------------------------------------------------------------
00086 snl_fei_tester::~snl_fei_tester()
00087 {
00088   delete linSysCore_;
00089   delete feData_;
00090 }
00091 
00092 //----------------------------------------------------------------------------
00093 int snl_fei_tester::testInitialization()
00094 {
00095   if (factory_.get() == NULL) {
00096     try {
00097       factory_ = fei::create_fei_Factory(comm_, data_->solverLibraryName_.c_str());
00098     }
00099     catch (std::runtime_error& exc) {
00100       fei::console_out() << exc.what()<<FEI_ENDL;
00101       return(1);
00102     }
00103     if (factory_.get() == NULL) {
00104       ERReturn(-1);
00105     }
00106   }
00107 
00108   std::vector<std::string> stdstrings;
00109   fei::utils::char_ptrs_to_strings(data_->numParams_, data_->paramStrings_,
00110                                    stdstrings);
00111   fei::ParameterSet paramset;
00112   fei::utils::parse_strings(stdstrings, " ", paramset);
00113 
00114   if (!path_.empty()) {
00115     paramset.add(fei::Param("debugOutput", path_.c_str()));
00116   }
00117 
00118   factory_->parameters(paramset);
00119 
00120   vecSpace_ = factory_->createVectorSpace(comm_, NULL);
00121 
00122   vecSpace_->setParameters(paramset);
00123 
00124   defineFieldsAndIDTypes();
00125 
00126   fei::SharedPtr<fei::VectorSpace> dummy;
00127   matrixGraph_ = factory_->createMatrixGraph(vecSpace_, dummy, NULL);
00128 
00129   matrixGraph_->setParameters(paramset);
00130 
00131   CHK_ERR( initElemBlocks() );
00132 
00133   CHK_ERR( initConstraints() );
00134 
00135   int i;
00136   for(i=0; i<data_->numSharedNodeSets_; ++i) {
00137     CommNodeSet& nodeSet = data_->sharedNodeSets_[i];
00138 
00139     CHK_ERR( vecSpace_->initSharedIDs(nodeSet.numNodes_,
00140                                        idTypes_[nodeTypeOffset_],
00141                                        nodeSet.nodeIDs_,
00142                                        nodeSet.procsPerNode_,
00143                                        nodeSet.procs_) );
00144   }
00145 
00146   CHK_ERR( matrixGraph_->initComplete() );
00147 
00148   return(0);
00149 }
00150 
00151 //----------------------------------------------------------------------------
00152 void snl_fei_tester::dumpMatrixFiles()
00153 {
00154   FEI_OSTRINGSTREAM osstr;
00155   osstr << "A_" << A_->typeName() << ".np"<<numProcs_;
00156   std::string str = osstr.str();
00157   A_->writeToFile(str.c_str());
00158 }
00159 
00160 //----------------------------------------------------------------------------
00161 void snl_fei_tester::setParameter(const char* param)
00162 {
00163   std::vector<std::string> stdstrings;
00164   fei::utils::char_ptrs_to_strings(1, &param, stdstrings);
00165   fei::ParameterSet paramset;
00166   fei::utils::parse_strings(stdstrings, " ", paramset);
00167   factory_->parameters(paramset);
00168   vecSpace_->setParameters(paramset);
00169   matrixGraph_->setParameters(paramset);
00170 
00171   linSys_->parameters(1, &param);
00172   A_->parameters(paramset);
00173 }
00174 
00175 //----------------------------------------------------------------------------
00176 int snl_fei_tester::testLoading()
00177 {
00178   linSys_ = factory_->createLinearSystem(matrixGraph_);
00179 
00180   A_ = factory_->createMatrix(matrixGraph_);
00181   x_ = factory_->createVector(matrixGraph_, true);
00182   b_ = factory_->createVector(matrixGraph_);
00183 
00184   matrixGraph_->setIndicesMode(fei::MatrixGraph::POINT_ENTRY_GRAPH);
00185 
00186   CHK_ERR( linSys_->parameters(data_->numParams_, data_->paramStrings_) );
00187 
00188   std::vector<std::string> stdstrings;
00189   fei::utils::char_ptrs_to_strings(data_->numParams_, data_->paramStrings_, stdstrings);
00190   fei::ParameterSet paramset;
00191   fei::utils::parse_strings(stdstrings, " ", paramset);
00192   CHK_ERR( A_->parameters(paramset) );
00193 
00194   linSys_->setMatrix(A_);
00195   linSys_->setRHS(b_);
00196   linSys_->setSolutionVector(x_);
00197 
00198   CHK_ERR( A_->putScalar(0.0) );
00199   CHK_ERR( b_->putScalar(0.0) );
00200 
00201   matrixGraph_->createSlaveMatrices();
00202 
00203   CHK_ERR( loadElemBlocks() );
00204 
00205   CHK_ERR( loadConstraints() );
00206 
00207   int i;
00208   for(i=0; i<data_->numBCNodeSets_; ++i) {
00209     BCNodeSet& bcSet = data_->bcNodeSets_[i];
00210     int fieldSize = data_->getFieldSize(bcSet.fieldID_);
00211     if (fieldSize < 1) {
00212       continue;
00213     }
00214 
00215     CHK_ERR( linSys_->loadEssentialBCs(bcSet.numNodes_,
00216                                        bcSet.nodeIDs_,
00217                                        idTypes_[nodeTypeOffset_],
00218                                        bcSet.fieldID_,
00219                                        bcSet.offsetsIntoField_,
00220                                        bcSet.prescribed_values_) );
00221   }
00222 
00223   CHK_ERR( linSys_->loadComplete() );
00224 
00225   return(0);
00226 }
00227 
00228 //----------------------------------------------------------------------------
00229 int snl_fei_tester::testSolve()
00230 {
00231   fei::SharedPtr<fei::Solver> solver = factory_->createSolver();
00232 
00233   std::vector<std::string> stdstrings;
00234   fei::utils::char_ptrs_to_strings(data_->numParams_, data_->paramStrings_, stdstrings);
00235   fei::ParameterSet paramset;
00236   fei::utils::parse_strings(stdstrings," ",paramset);
00237 
00238   int status, itersTaken = 0;
00239   CHK_ERR( solver->solve(linSys_.get(),
00240                          NULL, //preconditioningMatrix
00241                          paramset, itersTaken, status) );
00242 
00243   CHK_ERR( x_->scatterToOverlap() );
00244 
00245   return(0);
00246 }
00247 
00248 //----------------------------------------------------------------------------
00249 int snl_fei_tester::testCheckResult()
00250 {
00251   CHK_ERR( save_block_node_soln(*data_, x_.get(), data_->solnFileName_.c_str(),
00252                                 numProcs_, localProc_, 1));
00253 
00254   CHK_ERR( save_block_elem_soln(*data_, x_.get(), data_->solnFileName_.c_str(),
00255                                 numProcs_, localProc_, 1));
00256 
00257   CHK_ERR( save_multiplier_soln(*data_, x_.get(), data_->solnFileName_.c_str(),
00258                                 numProcs_, localProc_, 1));
00259 
00260   int err = SolnCheck::checkSolution(localProc_, numProcs_, data_->solnFileName_.c_str(),
00261                                      data_->checkFileName_.c_str(), "node", 1);
00262 
00263   err += SolnCheck::checkSolution(localProc_, numProcs_, data_->solnFileName_.c_str(),
00264                                   data_->checkFileName_.c_str(), "elem", 1);
00265 
00266   err += SolnCheck::checkSolution(localProc_, numProcs_, data_->solnFileName_.c_str(),
00267                                   data_->checkFileName_.c_str(), "mult", 1);
00268   int globalErr = err;
00269 #ifndef FEI_SER
00270   if (MPI_SUCCESS != MPI_Allreduce(&err, &globalErr, 1, MPI_INT, MPI_SUM,
00271                                    comm_)) return(-1);
00272 #endif
00273   if (globalErr != 0) return(-1);
00274   return(0);
00275 }
00276 
00277 //----------------------------------------------------------------------------
00278 void snl_fei_tester::defineFieldsAndIDTypes()
00279 {
00280   vecSpace_->defineFields(data_->numFields_, data_->fieldIDs_, data_->fieldSizes_);
00281 
00282   //nodeIDType == 0
00283   idTypes_.push_back(0);
00284 
00285   //constraintIDType == 1
00286   idTypes_.push_back(1);
00287 
00288   //elemDofIDType == 2
00289   idTypes_.push_back(2);
00290 
00291   vecSpace_->defineIDTypes(idTypes_.size(), &idTypes_[0] );
00292 
00293   nodeTypeOffset_ = 0;
00294   constraintTypeOffset_ = 1;
00295   elemTypeOffset_ = 2;
00296 }
00297 
00298 //----------------------------------------------------------------------------
00299 int snl_fei_tester::initElemBlocks()
00300 {
00301   for(int i=0; i<data_->numElemBlocks_; ++i) {
00302     ElemBlock& eb = data_->elemBlocks_[i];
00303 
00304     int patternID;
00305     definePattern(eb, patternID);
00306 
00307     CHK_ERR( matrixGraph_->initConnectivityBlock(eb.blockID_,
00308                                                eb.numElements_,
00309                                                patternID) );
00310 
00311     for(int j=0; j<eb.numElements_; ++j) {
00312       std::vector<int> conn(eb.numNodesPerElement_);
00313       for(int ii=0; ii<eb.numNodesPerElement_; ++ii) {
00314         conn[ii] = eb.elemConn_[j][ii];
00315       }
00316 
00317       CHK_ERR( matrixGraph_->initConnectivity(eb.blockID_,
00318                                               eb.elemIDs_[j],
00319                                               &conn[0]) );
00320     }
00321   }
00322 
00323   return(0);
00324 }
00325 
00326 //----------------------------------------------------------------------------
00327 int snl_fei_tester::loadElemBlocks()
00328 {
00329   int i;
00330   for(i=0; i<data_->numElemBlocks_; ++i) {
00331     ElemBlock& eb = data_->elemBlocks_[i];
00332 
00333     if (eb.numElements_ < 1) {
00334       continue;
00335     }
00336 
00337     int numIndices = matrixGraph_->getConnectivityNumIndices(eb.blockID_);
00338 
00339     std::vector<int> indices(numIndices);
00340 
00341     for(int j=0; j<eb.numElements_; ++j) {
00342       int checkNum;
00343       CHK_ERR( matrixGraph_->getConnectivityIndices(eb.blockID_,
00344                                                     eb.elemIDs_[j],
00345                                                     numIndices,
00346                                                     &indices[0],
00347                                                     checkNum) );
00348       if (numIndices != checkNum) {
00349         ERReturn(-1);
00350       }
00351 
00352       CHK_ERR( A_->sumIn(eb.blockID_, eb.elemIDs_[j],
00353                          eb.elemStiff_[j]) );
00354 
00355       CHK_ERR( b_->sumIn(numIndices, &indices[0],
00356                          eb.elemLoad_[j], 0) );
00357     }
00358   }
00359 
00360   return(0);
00361 }
00362 
00363 //----------------------------------------------------------------------------
00364 int snl_fei_tester::initConstraints()
00365 {
00366   std::vector<int> idTypes;
00367   int constraintID = localProc_*100000;
00368   int i;
00369   for(i=0; i<data_->numCRMultSets_; ++i) {
00370     CRSet& crSet = data_->crMultSets_[i];
00371 
00372     for(int j=0; j<1; ++j) {
00373       idTypes.assign(crSet.numNodes_, idTypes_[nodeTypeOffset_]);
00374 
00375       crSet.crID_ = constraintID++;
00376       int constraintIDType = idTypes_[constraintTypeOffset_];
00377       CHK_ERR( matrixGraph_->initLagrangeConstraint(crSet.crID_,
00378                                                   constraintIDType,
00379                                                   crSet.numNodes_,
00380                                                   &idTypes[0],
00381                                                   crSet.nodeIDs_[j],
00382                                                   crSet.fieldIDs_) );
00383     }
00384   }
00385 
00386   for(i=0; i<data_->numCRPenSets_; ++i) {
00387     CRSet& crSet = data_->crPenSets_[i];
00388 
00389     for(int j=0; j<1; ++j) {
00390       idTypes.assign(crSet.numNodes_, idTypes_[nodeTypeOffset_]);
00391 
00392       crSet.crID_ = constraintID++;
00393       int constraintIDType = idTypes_[constraintTypeOffset_];
00394       CHK_ERR( matrixGraph_->initPenaltyConstraint(crSet.crID_,
00395                                                   constraintIDType,
00396                                                   crSet.numNodes_,
00397                                                   &idTypes[0],
00398                                                   crSet.nodeIDs_[j],
00399                                                   crSet.fieldIDs_) );
00400     }
00401   }
00402 
00403   std::map<int,int> fieldDB;
00404   for(i=0; i<data_->numFields_; ++i) {
00405     fieldDB.insert(std::pair<int,int>(data_->fieldIDs_[i], data_->fieldSizes_[i]));
00406   }
00407 
00408   std::vector<int> nodeIDs;
00409   std::vector<int> fieldIDs;
00410   std::vector<double> weights;
00411 
00412   for(i=0; i<data_->numSlaveVars_; i++) {
00413     int ii;
00414     CRSet& crSet = data_->slaveVars_[i];
00415 
00416     nodeIDs.resize(crSet.numNodes_+1);
00417     nodeIDs[0] = crSet.slaveNodeID_;
00418     fieldIDs.resize(0);
00419     fieldIDs.push_back(crSet.slaveFieldID_);
00420 
00421     for(ii=0; ii<crSet.numNodes_; ++ii) {
00422       nodeIDs[ii+1] = crSet.nodeIDs_[0][ii];
00423       fieldIDs.push_back(crSet.fieldIDs_[ii]);
00424     }
00425 
00426     idTypes.assign(crSet.numNodes_+1, idTypes_[nodeTypeOffset_]);
00427 
00428     int fieldSize = fieldDB[crSet.slaveFieldID_];
00429     weights.resize(0);
00430     for(ii=0; ii<fieldSize; ++ii) weights.push_back(0.0);
00431     weights[crSet.slaveOffset_] = -1.0;
00432     int offset = 0;
00433     for(ii=0; ii<crSet.numNodes_; ++ii) {
00434       fieldSize = fieldDB[crSet.fieldIDs_[ii]];
00435       for(int jj=0; jj<fieldSize; ++jj) {
00436         weights.push_back(crSet.weights_[offset++]);
00437       }
00438     }
00439 
00440     CHK_ERR( matrixGraph_->initSlaveConstraint(crSet.numNodes_+1,
00441                                                &idTypes[0],
00442                                                &nodeIDs[0],
00443                                                &fieldIDs[0],
00444                                                0,
00445                                                crSet.slaveOffset_,
00446                                                &weights[0],
00447                                                crSet.values_[0]));
00448   }
00449 
00450   return(0);
00451 }
00452 
00453 //----------------------------------------------------------------------------
00454 int snl_fei_tester::loadConstraints()
00455 {
00456   int i;
00457   for(i=0; i<data_->numCRMultSets_; ++i) {
00458     CRSet& crSet = data_->crMultSets_[i];
00459 
00460     for(int j=0; j<1; ++j) {
00461       CHK_ERR( linSys_->loadLagrangeConstraint(crSet.crID_,
00462                                                crSet.weights_,
00463                                                crSet.values_[j]) );
00464     }
00465   }
00466 
00467   for(i=0; i<data_->numCRPenSets_; ++i) {
00468     CRSet& crSet = data_->crPenSets_[i];
00469 
00470     for(int j=0; j<1; ++j) {
00471       CHK_ERR( linSys_->loadPenaltyConstraint(crSet.crID_,
00472                                               crSet.weights_,
00473                                               crSet.penValues_[j],
00474                                               crSet.values_[j]) );
00475     }
00476   }
00477 
00478   return(0);
00479 }
00480 
00481 //----------------------------------------------------------------------------
00482 void snl_fei_tester::definePattern(ElemBlock& eb, int& patternID)
00483 {
00484   int i, j, numIDTypes = 1;
00485   numIDTypes += eb.numElemDOF_>0 ? 1 : 0;
00486 
00487   //find out how many nodal fields there are, total.
00488   std::vector<int> nodalFieldIDs;
00489   std::vector<int> flatFieldIDsArray;
00490   for(i=0; i<eb.numNodesPerElement_; ++i) {
00491     for(j=0; j<eb.numFieldsPerNode_[i]; ++j) {
00492       fei::sortedListInsert(eb.nodalFieldIDs_[i][j], nodalFieldIDs);
00493       flatFieldIDsArray.push_back(eb.nodalFieldIDs_[i][j]);
00494     }
00495   }
00496 
00497   patternID = numPatterns_++;
00498 
00499   if (numIDTypes == 1 && nodalFieldIDs.size() == 1) {
00500     //This is a very simple pattern
00501     patternID = matrixGraph_->definePattern(eb.numNodesPerElement_,
00502                                 idTypes_[nodeTypeOffset_],
00503                                 nodalFieldIDs[0]);
00504   }
00505   else if (numIDTypes == 1) {
00506     std::vector<int> numFieldsPerID(eb.numNodesPerElement_);
00507 
00508     patternID = matrixGraph_->definePattern(eb.numNodesPerElement_,
00509                                 idTypes_[nodeTypeOffset_],
00510                                 eb.numFieldsPerNode_,
00511                                 &flatFieldIDsArray[0]);
00512   }
00513   else {
00514     std::vector<int> idTypes(eb.numNodesPerElement_+1, idTypes_[nodeTypeOffset_]);
00515     idTypes[idTypes.size()-1] = idTypes_[elemTypeOffset_];
00516     std::vector<int> numFieldsPerID(idTypes.size());
00517     std::vector<int> fieldIDs;
00518     for(i=0; i<eb.numNodesPerElement_; ++i) {
00519       numFieldsPerID[i] = eb.numFieldsPerNode_[i];
00520       for(j=0; j<eb.numFieldsPerNode_[i]; ++j) {
00521         fieldIDs.push_back(eb.nodalFieldIDs_[i][j]);
00522       }
00523     }
00524     numFieldsPerID[idTypes.size()-1] = eb.numElemDOF_;
00525     for(i=0; i<eb.numElemDOF_; ++i) {
00526       fieldIDs.push_back(eb.elemDOFFieldIDs_[i]);
00527     }
00528 
00529     patternID = matrixGraph_->definePattern(idTypes.size(),
00530                                 &idTypes[0],
00531                                 &numFieldsPerID[0],
00532                                 &fieldIDs[0]);
00533   }
00534 }
00535 
00536 //----------------------------------------------------------------------------
00537 int snl_fei_tester::save_block_node_soln(DataReader& data, fei::Vector* vec,
00538                                      const char* solnFileName, int numProcs,
00539                                      int localProc, int solveCounter)
00540 {
00541   (void)solveCounter;
00542 
00543   int numLocalNodes = vecSpace_->getNumOwnedAndSharedIDs(idTypes_[nodeTypeOffset_]);
00544 
00545   int* nodeList = new int[numLocalNodes];
00546 
00547   int checkNum = 0;
00548   int err = vecSpace_->getOwnedAndSharedIDs(idTypes_[nodeTypeOffset_],
00549                                     numLocalNodes, nodeList, checkNum);
00550   if (err != 0) {
00551     ERReturn(-1);
00552   }
00553 
00554   FEI_OSTRINGSTREAM fileName;
00555   fileName<< solnFileName<<".node."<<solveCounter<<"."<<numProcs<<"."<<localProc;
00556   std::string str = fileName.str();
00557   FEI_OFSTREAM outfile(str.c_str());
00558 
00559   if (!outfile || outfile.bad()) {
00560     fei::console_out() << "ERROR opening solution output file " << fileName << FEI_ENDL;
00561     return(-1);
00562   }
00563 
00564   outfile.setf(IOS_SCIENTIFIC, IOS_FLOATFIELD);
00565 
00566   std::vector<double> solnData;
00567   std::vector<int> fieldList;
00568 
00569   int totalSize = 0;
00570 
00571   for(int i=0; i<numLocalNodes; i++) {
00572     int idType = idTypes_[nodeTypeOffset_];
00573     int ID = nodeList[i];
00574 
00575     int numDOF = vecSpace_->getNumDegreesOfFreedom(idType, ID);
00576     int numFields = vecSpace_->getNumFields(idType, ID);
00577     solnData.resize(numDOF);
00578     vecSpace_->getFields(idType, ID, fieldList);
00579 
00580     outfile << ID << " " << numDOF << FEI_ENDL;
00581     for(int j=0; j<numFields; ++j) {
00582       int fieldSize = vecSpace_->getFieldSize(fieldList[j]);
00583       totalSize += fieldSize;
00584 
00585       CHK_ERR( vec->copyOutFieldData(fieldList[j], idType,
00586                                      1, &ID, &solnData[0]) );
00587 
00588       for(int k=0; k<fieldSize; ++k) {
00589         outfile << solnData[k] << " ";
00590       }
00591     }
00592     outfile << FEI_ENDL;
00593   }
00594 
00595   FEI_COUT << "save-node-soln: wrote " << totalSize << " entries for " << numLocalNodes << " nodes to " << str << FEI_ENDL;
00596 
00597   delete [] nodeList;
00598 
00599   outfile.close();
00600   return(0);
00601 }
00602 
00603 //----------------------------------------------------------------------------
00604 int snl_fei_tester::save_block_elem_soln(DataReader& data, fei::Vector* vec,
00605                                          const char* solnFileName,
00606                                          int numProcs,
00607                                          int localProc, int solveCounter)
00608 {
00609   (void)solveCounter;
00610 
00611   int numLocalElems = vecSpace_->getNumOwnedAndSharedIDs(idTypes_[elemTypeOffset_]);
00612 
00613   int* elemList = new int[numLocalElems];
00614 
00615   int checkNum = 0;
00616   int err = vecSpace_->getOwnedAndSharedIDs(idTypes_[elemTypeOffset_],
00617                                     numLocalElems, elemList, checkNum);
00618   if (err != 0) {
00619     ERReturn(-1);
00620   }
00621 
00622   FEI_OSTRINGSTREAM fileName;
00623   fileName<< solnFileName<<".elem."<<solveCounter<<"."<<numProcs<<"."<<localProc;
00624   std::string str = fileName.str();
00625   FEI_OFSTREAM outfile(str.c_str());
00626 
00627   if (!outfile || outfile.bad()) {
00628     fei::console_out() << "ERROR opening solution output file " << fileName << FEI_ENDL;
00629     return(-1);
00630   }
00631 
00632   std::vector<double> solnData;
00633   std::vector<int> fieldList;
00634 
00635   for(int i=0; i<numLocalElems; i++) {
00636     int idType = idTypes_[elemTypeOffset_];
00637     int ID = elemList[i];
00638 
00639     int numDOF = vecSpace_->getNumDegreesOfFreedom(idType, ID);
00640     int numFields = vecSpace_->getNumFields(idType, ID);
00641     solnData.resize(numDOF);
00642     vecSpace_->getFields(idType, ID, fieldList);
00643 
00644     outfile << ID << " " << numDOF << FEI_ENDL;
00645     for(int j=0; j<numFields; ++j) {
00646       int fieldSize = vecSpace_->getFieldSize(fieldList[j]);
00647 
00648       CHK_ERR( vec->copyOutFieldData(fieldList[j], idType,
00649                                      1, &ID, &solnData[0]) );
00650 
00651       for(int k=0; k<fieldSize; ++k) {
00652         outfile << solnData[k] << " ";
00653       }
00654     }
00655     outfile << FEI_ENDL;
00656   }
00657 
00658   delete [] elemList;
00659 
00660   outfile.close();
00661   return(0);
00662 }
00663 
00664 //----------------------------------------------------------------------------
00665 int snl_fei_tester::save_multiplier_soln(DataReader& data, fei::Vector* vec,
00666                                          const char* solnFileName,
00667                                          int numProcs, int localProc,
00668                                          int solveCounter)
00669 {
00670   (void)solveCounter;
00671 
00672   int numLocalCRs = vecSpace_->getNumOwnedAndSharedIDs(idTypes_[constraintTypeOffset_]);
00673 
00674    int* globalNumCRs = new int[numProcs];
00675 #ifndef FEI_SER
00676    if (MPI_Allgather(&numLocalCRs, 1, MPI_INT, globalNumCRs, 1, MPI_INT,
00677                  comm_) != MPI_SUCCESS) {
00678      ERReturn(-1);
00679    }
00680 #endif
00681 
00682    int localCRStart = 0;
00683 #ifndef FEI_SER
00684    for(int p=0; p<localProc; p++) localCRStart += globalNumCRs[p];
00685 #endif
00686 
00687    delete [] globalNumCRs;
00688 
00689   std::vector<int> crList(numLocalCRs);
00690 
00691   int checkNum = 0;
00692   int err = vecSpace_->getOwnedAndSharedIDs(
00693     idTypes_[constraintTypeOffset_], numLocalCRs,
00694     numLocalCRs ? &crList[0] : 0, checkNum);
00695   if (err != 0) {
00696     ERReturn(-1);
00697   }
00698 
00699   FEI_OSTRINGSTREAM fileName;
00700   fileName<< solnFileName<<".mult."<<solveCounter<<"."<<numProcs<<"."<<localProc;
00701   std::string str = fileName.str();
00702   FEI_OFSTREAM outfile(str.c_str());
00703 
00704   if (!outfile || outfile.bad()) {
00705     fei::console_out() << "ERROR opening solution output file " << fileName << FEI_ENDL;
00706     return(-1);
00707   }
00708 
00709   std::vector<double> solnData;
00710   std::vector<int> fieldList;
00711 
00712   for(int i=0; i<numLocalCRs; i++) {
00713     int idType = idTypes_[constraintTypeOffset_];
00714     int ID = crList[i];
00715 
00716     solnData.resize(1);
00717 
00718     outfile << localCRStart++ << " " << 1 << FEI_ENDL;
00719     for(int j=0; j<1; ++j) {
00720       int globalIndex = -1;
00721       CHK_ERR( vecSpace_->getGlobalIndex(idType, ID, globalIndex) );
00722 
00723       CHK_ERR( vec->copyOut(1, &globalIndex, &solnData[0]) );
00724 
00725       for(int k=0; k<1; ++k) {
00726         outfile << solnData[k] << " ";
00727       }
00728     }
00729     outfile << FEI_ENDL;
00730   }
00731 
00732   outfile.close();
00733   return(0);
00734 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends