FEI_tester.cpp

00001 /*--------------------------------------------------------------------*/
00002 /*    Copyright 2005 Sandia Corporation.                              */
00003 /*    Under the terms of Contract DE-AC04-94AL85000, there is a       */
00004 /*    non-exclusive license for use of this work by or on behalf      */
00005 /*    of the U.S. Government.  Export of this program may require     */
00006 /*    a license from the United States Government.                    */
00007 /*--------------------------------------------------------------------*/
00008 
00009 #include <fei_sstream.hpp>
00010 #include <fei_fstream.hpp>
00011 
00012 #include <test_utils/fei_test_utils.hpp>
00013 
00014 #include <test_utils/FEI_tester.hpp>
00015 
00016 #include <fei_LinearSystemCore.hpp>
00017 #include <fei_LibraryWrapper.hpp>
00018 #include <snl_fei_Utils.hpp>
00019 
00020 #include <fei_FEI_Impl.hpp>
00021 
00022 #include <test_utils/LibraryFactory.hpp>
00023 
00024 #ifdef HAVE_FEI_FETI
00025 #include <FETI_DP_FiniteElementData.h>
00026 #endif
00027 
00028 #include <test_utils/DataReader.hpp>
00029 #include <test_utils/SolnCheck.hpp>
00030 
00031 #undef fei_file
00032 #define fei_file "FEI_tester.cpp"
00033 
00034 #include <fei_ErrMacros.hpp>
00035 
00036 //----------------------------------------------------------------------------
00037 FEI_tester::FEI_tester(fei::SharedPtr<DataReader> data_reader,
00038            MPI_Comm comm, int localProc, int numProcs, bool useNewFEI)
00039   : comm_(comm),
00040     fei_(),
00041     wrapper_(),
00042     data_(data_reader),
00043     localProc_(localProc),
00044     numProcs_(numProcs),
00045     numMatrices(1),
00046     matrixIDs(NULL),
00047     numRHSs(1),
00048     rhsIDs(NULL),
00049     useNewFEI_(useNewFEI)
00050 {
00051 }
00052 
00053 //----------------------------------------------------------------------------
00054 FEI_tester::~FEI_tester()
00055 {
00056   delete [] matrixIDs;
00057   delete [] rhsIDs;
00058 }
00059 
00060 //----------------------------------------------------------------------------
00061 int FEI_tester::testInitialization()
00062 {
00063   if (data_.get() == NULL) {
00064     ERReturn(-1);
00065   }
00066 
00067   CHK_ERR( createFEIinstance(data_->solverLibraryName_.c_str()) );
00068 
00069   const char* feiVersionString;
00070   CHK_ERR( fei_->version(feiVersionString) );
00071 
00072   FEI_COUT << "FEI version: " << feiVersionString << FEI_ENDL;
00073 
00074   fei_->parameters(data_->numParams_, data_->paramStrings_);
00075 
00076   //Now we check the solveType. A regular Ax=b solve corresponds to 
00077   //solveType==FEI_SINGLE_SYSTEM. The aggregate stuff (solveType==
00078   //FEI_AGGREGATE_SUM) is for power users who
00079   //want to assemble more than one matrix system and solve a linear
00080   //combination of them, an aggregate system.
00081 
00082   if (data_->solveType_ == FEI_AGGREGATE_SUM) {
00083     CHK_ERR( setIDlists());
00084   }
00085 
00086   CHK_ERR( initializationPhase() );
00087 
00088   int numBlkActNodes;
00089   for(int i=0; i<data_->numElemBlocks_; ++i) {
00090     ElemBlock& eblk = data_->elemBlocks_[i];
00091     int elemBlockID = eblk.blockID_;
00092     CHK_ERR( fei_->getNumBlockActNodes(elemBlockID, numBlkActNodes) );
00093   }
00094 
00095   //************** Initialization Phase is now complete *****************
00096 
00097   return(0);
00098 }
00099 
00100 //----------------------------------------------------------------------------
00101 int FEI_tester::testLoading()
00102 {
00103   CHK_ERR(fei_->resetRHSVector());
00104   CHK_ERR(fei_->resetMatrix());    //try out all of these reset functions,
00105   CHK_ERR(fei_->resetSystem());    //just for the sake of coverage...
00106 
00107   // ************ Load Phase (delegated to a function) ***************
00108 
00109   //obviously only one of these load-phases should be
00110   //performed.
00111 
00112   if (data_->solveType_ == FEI_SINGLE_SYSTEM) {
00113     CHK_ERR( normalLoadPhase());
00114   }
00115   if (data_->solveType_ == FEI_AGGREGATE_SUM) {
00116     CHK_ERR( aggregateLoadPhase());
00117   }
00118 
00119   CHK_ERR( fei_->loadComplete() );
00120 
00121   CHK_ERR( exerciseResidualNorm() );
00122 
00123   //**** let's try out the 'put' functions. ******************
00124   CHK_ERR( exercisePutFunctions() );
00125 
00126   return(0);
00127 }
00128 
00129 //----------------------------------------------------------------------------
00130 int FEI_tester::testSolve()
00131 {
00132   //If solverLibraryName is TEST_LSC, then we're not running a real solver so
00133   //just return.
00134   std::string sname(data_->solverLibraryName_);
00135   if (sname == "TEST_LSC") {
00136     return( 0 );
00137   }
00138 
00139   int status;
00140   int err = fei_->solve(status);
00141 
00142   //FEI_COUT << "fei->solve, err = " << err << ", status = " << status << FEI_ENDL;
00143 
00144   if (err != 0 || status != 0) {
00145     FEI_COUT << "!!!! solve returned: err: "<<err<<", status: "<<status<<FEI_ENDL;
00146     return(err);
00147   }
00148 
00149   if (localProc_ == 0) {
00150     int itersTaken;
00151     CHK_ERR( fei_->iterations(itersTaken));
00152     //FEI_COUT << "iterations: " << itersTaken << FEI_ENDL;
00153   }
00154 
00155   CHK_ERR( exerciseResidualNorm() );
00156 
00157   return(0);
00158 }
00159 
00160 //----------------------------------------------------------------------------
00161 void FEI_tester::dumpMatrixFiles()
00162 {
00163 }
00164 
00165 //----------------------------------------------------------------------------
00166 void FEI_tester::setParameter(const char*)
00167 {
00168 }
00169 
00170 //----------------------------------------------------------------------------
00171 int FEI_tester::testCheckResult()
00172 {
00173   //If solverLibraryName is TEST_LSC, then we're not running a real solver so
00174   //just check the matrix.
00175   std::string sname(data_->solverLibraryName_);
00176   if (sname == "TEST_LSC") {
00177     return( lsc_matrix_check() );
00178   }
00179 
00180   CHK_ERR( save_block_node_soln(*data_, *fei_, data_->solnFileName_.c_str(),
00181         numProcs_, localProc_, 1));
00182 
00183   CHK_ERR( save_block_elem_soln(*data_, *fei_, data_->solnFileName_.c_str(),
00184         numProcs_, localProc_, 1));
00185 
00186   CHK_ERR( save_multiplier_soln(*data_, *fei_, data_->solnFileName_.c_str(),
00187         numProcs_, localProc_, 1));
00188 
00189   int err = SolnCheck::checkSolution(localProc_, numProcs_, data_->solnFileName_.c_str(),
00190              data_->checkFileName_.c_str(), "node", 1);
00191 
00192   err += SolnCheck::checkSolution(localProc_, numProcs_, data_->solnFileName_.c_str(),
00193           data_->checkFileName_.c_str(), "elem", 1);
00194 
00195   err += SolnCheck::checkSolution(localProc_, numProcs_, data_->solnFileName_.c_str(),
00196           data_->checkFileName_.c_str(), "mult", 1);
00197   int globalErr = err;
00198 #ifndef FEI_SER
00199   if (MPI_SUCCESS != MPI_Allreduce(&err, &globalErr, 1, MPI_INT, MPI_SUM,
00200            comm_)) return(-1);
00201 #endif
00202   if (globalErr != 0) {
00203     ERReturn(-1);
00204   }
00205 
00206   return(0);
00207 }
00208 
00209 //----------------------------------------------------------------------------
00210 int FEI_tester::createFEIinstance(const char* solverName)
00211 {
00212   try {
00213     wrapper_ = fei::create_LibraryWrapper(comm_, solverName);
00214   }
00215   catch(std::runtime_error& exc) {
00216     FEI_CERR << exc.what()<<FEI_ENDL;
00217     return(1);
00218   }
00219 
00220   if (wrapper_.get() == NULL) ERReturn(-1);
00221 
00222   if (useNewFEI_) {
00223     fei_.reset(new fei::FEI_Impl(wrapper_, comm_, 0));
00224   }
00225   else {
00226     fei_.reset(new FEI_Implementation(wrapper_, comm_, 0));
00227   }
00228 
00229   if (fei_.get() == NULL) ERReturn(-1);
00230 
00231   return(0);
00232 }
00233 
00234 //----------------------------------------------------------------------------
00235 int FEI_tester::setIDlists()
00236 {
00237   snl_fei::getIntParamValue("numMatrices",
00238           data_->numParams_,
00239           data_->paramStrings_,
00240           numMatrices);
00241 
00242   matrixIDs = new int[numMatrices];
00243   numRHSs = 1;
00244   rhsIDs = new int[1];
00245   rhsIDs[0] = 0;
00246 
00247   for(int i=0; i<numMatrices; i++) {
00248     matrixIDs[i] = i;
00249   }
00250 
00251   CHK_ERR(fei_->setIDLists(numMatrices, matrixIDs, numRHSs, rhsIDs));
00252   return(0);
00253 }
00254 
00255 //----------------------------------------------------------------------------
00256 int FEI_tester::initializationPhase()
00257 {
00258   if (data_->solveType_ != FEI_SINGLE_SYSTEM &&
00259       data_->solveType_ != FEI_AGGREGATE_SUM) {
00260     FEI_COUT << "FEI_tester: bad solveType: " << data_->solveType_ << FEI_ENDL;
00261     return(-1);
00262   }
00263 
00264   CHK_ERR( fei_->setSolveType(data_->solveType_));
00265 
00266   CHK_ERR(fei_->initFields(data_->numFields_, data_->fieldSizes_, data_->fieldIDs_));
00267 
00268   int i;
00269   for(i=0; i<data_->numElemBlocks_; i++) {
00270     ElemBlock& block = data_->elemBlocks_[i];
00271 
00272     CHK_ERR(fei_->initElemBlock(block.blockID_,
00273              block.numElements_,
00274              block.numNodesPerElement_,
00275              block.numFieldsPerNode_,
00276              block.nodalFieldIDs_,
00277              block.numElemDOF_,
00278              block.elemDOFFieldIDs_,
00279              block.interleaveStrategy_) );
00280 
00281     for(int el=0; el<block.numElements_; el++) {
00282       CHK_ERR(fei_->initElem(block.blockID_, 
00283           block.elemIDs_[el],
00284           block.elemConn_[el]));
00285     }
00286   }
00287 
00288   for(i=0; i<data_->numSharedNodeSets_; i++) {
00289     CommNodeSet& shNodeSet = data_->sharedNodeSets_[i];
00290 
00291     CHK_ERR(fei_->initSharedNodes(shNodeSet.numNodes_, shNodeSet.nodeIDs_,
00292          shNodeSet.procsPerNode_, shNodeSet.procs_));
00293   }
00294 
00295   //********* Initialize any slave variables *****************************
00296 
00297   for(i=0; i<data_->numSlaveVars_; i++) {
00298     CRSet& crSet = data_->slaveVars_[i];
00299 
00300     CHK_ERR( fei_->initSlaveVariable(crSet.slaveNodeID_,
00301             crSet.slaveFieldID_,
00302             crSet.slaveOffset_,
00303             crSet.numNodes_,
00304             crSet.nodeIDs_[0],
00305             crSet.fieldIDs_,
00306             crSet.weights_,
00307             crSet.values_[0]) );
00308   }
00309 
00310   //*********** Initialize Constraint Relation Equations *****************
00311 
00312   for(i=0; i<data_->numCRMultSets_; i++) {
00313     CRSet& crSet = data_->crMultSets_[i];
00314 
00315     for(int j=0; j<1; j++) {
00316       CHK_ERR(fei_->initCRMult(crSet.numNodes_,
00317             crSet.nodeIDs_[j],
00318             crSet.fieldIDs_,
00319             crSet.crID_));
00320     }
00321   }
00322 
00323   for(i=0; i<data_->numCRPenSets_; i++) {
00324     CRSet& crSet = data_->crPenSets_[i];
00325 
00326     for(int j=0; j<1; j++) {
00327       CHK_ERR(fei_->initCRPen(crSet.numNodes_,
00328            crSet.nodeIDs_[j],
00329            crSet.fieldIDs_,
00330            crSet.crID_));
00331     }
00332   }
00333 
00334   CHK_ERR(fei_->initComplete());
00335   return(0);
00336 }
00337 
00338 //----------------------------------------------------------------------------
00339 int FEI_tester::normalLoadPhase()
00340 {
00341   int i;
00342 
00343   //*************** Boundary Condition Node Sets *************************
00344 
00345   for(i=0; i<data_->numBCNodeSets_; i++) {
00346     BCNodeSet& bcSet = data_->bcNodeSets_[i];
00347 
00348     CHK_ERR(fei_->loadNodeBCs(bcSet.numNodes_,
00349            bcSet.nodeIDs_,
00350            bcSet.fieldID_,
00351            bcSet.offsetsIntoField_,
00352            bcSet.prescribed_values_));
00353   }
00354 
00355    for(i=0; i<data_->numElemBlocks_; i++) {
00356       ElemBlock& block = data_->elemBlocks_[i];
00357 
00358       for(int el=0; el<block.numElements_; el++) {
00359 
00360          CHK_ERR(fei_->sumInElemMatrix(block.blockID_,
00361                                 block.elemIDs_[el],
00362                                 block.elemConn_[el],
00363                                 block.elemStiff_[el],
00364                                 block.elemFormat_));
00365       }
00366    }
00367 
00368    for(i=0; i<data_->numElemBlocks_; i++) {
00369       ElemBlock& block = data_->elemBlocks_[i];
00370 
00371       for(int el=0; el<block.numElements_; el++) {
00372 
00373          CHK_ERR(fei_->sumInElemRHS(block.blockID_,
00374                                 block.elemIDs_[el],
00375                                 block.elemConn_[el],
00376                                 block.elemLoad_[el]));
00377       }
00378    }
00379 
00380    //******** Load Constraint Relation Equations ***********************
00381 
00382    for(i=0; i<data_->numCRMultSets_; i++) {
00383       CRSet& crSet = data_->crMultSets_[i];
00384 
00385       for(int j=0; j<1; j++) {
00386          CHK_ERR(fei_->loadCRMult(crSet.crID_,
00387                                  crSet.numNodes_,
00388                                  crSet.nodeIDs_[j],
00389                                  crSet.fieldIDs_,
00390                                  crSet.weights_,
00391                                  crSet.values_[j]))
00392       }
00393    }
00394 
00395    for(i=0; i<data_->numCRPenSets_; i++) {
00396       CRSet& crSet = data_->crPenSets_[i];
00397 
00398       for(int j=0; j<1; j++) {
00399          CHK_ERR(fei_->loadCRPen(crSet.crID_,
00400                                 crSet.numNodes_,
00401                                 crSet.nodeIDs_[j],
00402                                 crSet.fieldIDs_,
00403                                 crSet.weights_,
00404                                 crSet.values_[j],
00405                                 crSet.penValues_[j]))
00406       }
00407    }
00408   return(0);
00409 }
00410 
00411 //----------------------------------------------------------------------------
00412 int FEI_tester::aggregateLoadPhase()
00413 {
00414    int i;
00415 
00416    for(i=0; i<numMatrices; i++) {
00417       CHK_ERR(fei_->setCurrentMatrix(matrixIDs[i]))
00418 
00419       for(int j=0; j<data_->numElemBlocks_; j++) {
00420          ElemBlock& block = data_->elemBlocks_[j];
00421 
00422          for(int el=0; el<block.numElements_; el++) {
00423 
00424             CHK_ERR(fei_->sumInElemMatrix(block.blockID_,
00425                                    block.elemIDs_[el],
00426                                    block.elemConn_[el],
00427                                    block.elemStiff_[el],
00428                                    block.elemFormat_))
00429          }
00430       }
00431    }
00432 
00433    for(i=0; i<numRHSs; i++) {
00434       CHK_ERR(fei_->setCurrentRHS(rhsIDs[i]))
00435 
00436       for(int j=0; j<data_->numElemBlocks_; j++) {
00437          ElemBlock& block = data_->elemBlocks_[j];
00438 
00439          for(int el=0; el<block.numElements_; el++) {
00440             CHK_ERR(fei_->sumInElemRHS(block.blockID_,
00441                                    block.elemIDs_[el],
00442                                    block.elemConn_[el],
00443                                    block.elemLoad_[el]))
00444          }
00445       }
00446    }
00447 
00448    //*************** Boundary Condition Node Sets *************************
00449 
00450    for(i=0; i<data_->numBCNodeSets_; i++) {
00451       BCNodeSet& bcSet = data_->bcNodeSets_[i];
00452 
00453       CHK_ERR(fei_->loadNodeBCs(bcSet.numNodes_,
00454                      bcSet.nodeIDs_,
00455                      bcSet.fieldID_,
00456                      bcSet.offsetsIntoField_,
00457                      bcSet.prescribed_values_))
00458    }
00459 
00460    double* matScalars = new double[numMatrices];
00461    for(i=0; i<numMatrices; i++) {
00462       matScalars[i] = 1.0;
00463    }
00464 
00465    int rhsScaleID = rhsIDs[0];
00466    double rhsScalar = 1.0;
00467 
00468    CHK_ERR(fei_->setMatScalars(numMatrices, matrixIDs, matScalars))
00469    CHK_ERR(fei_->setRHSScalars(1, &rhsScaleID, &rhsScalar))
00470 
00471    delete [] matScalars;
00472   return(0);
00473 }
00474 
00475 //----------------------------------------------------------------------------
00476 int FEI_tester::exerciseResidualNorm()
00477 {
00478   std::string sname(data_->solverLibraryName_);
00479   if (sname == "TEST_LSC") {
00480     return( 0 );
00481   }
00482   //FEI_COUT << "numFields: " << data_->numFields_ << FEI_ENDL;
00483   double* norms = new double[data_->numFields_];
00484   int *fields = new int[data_->numFields_];
00485   for(int i=0; i<data_->numFields_; ++i) {
00486     fields[i] = data_->fieldIDs_[i];
00487   }
00488 
00489   CHK_ERR( fei_->residualNorm(1, data_->numFields_, fields, norms) );
00490   //int n;
00491   //for(n=0; n<data_->numFields_; n++) {
00492   //  FEI_COUT << " field["<<n<<"]: " << fields[n]
00493   //    << ", 1-norm: " << norms[n] << FEI_ENDL;
00494   //}
00495 
00496   delete [] fields;
00497   delete [] norms;
00498   return(0);
00499 }
00500 
00501 //----------------------------------------------------------------------------
00502 int FEI_tester::exercisePutFunctions()
00503 {
00504    //let's exercise the putNodalFieldData function.
00505 
00506    int numNodes;
00507    CHK_ERR( fei_->getNumLocalNodes(numNodes) );
00508    std::vector<int> nodeIDs(numNodes);
00509    int* nodeIDsPtr = &nodeIDs[0];
00510    int checkNumNodes;
00511    CHK_ERR( fei_->getLocalNodeIDList(checkNumNodes, nodeIDsPtr, numNodes) );
00512 
00513    for(int i=0; i<data_->numFields_; ++i) {
00514      int fieldID = data_->fieldIDs_[i];
00515      int fieldSize = data_->fieldSizes_[i];
00516      std::vector<double> data(numNodes*fieldSize, 0.0001);
00517 
00518      CHK_ERR( fei_->putNodalFieldData(fieldID, numNodes, nodeIDsPtr,
00519               &data[0]) );
00520    }
00521 
00522    return(0);
00523 }
00524 
00525 //----------------------------------------------------------------------------
00526 int FEI_tester::save_block_node_soln(DataReader& data, FEI& fei,
00527              const char* solnFileName, int numProcs,
00528              int localProc, int solveCounter)
00529 {
00530   (void)solveCounter;
00531    int i;
00532 
00533    int maxNumEqnsPerNode = 0;
00534    for(i=0; i<data.numFields_; ++i) {
00535      maxNumEqnsPerNode += data.fieldSizes_[i];
00536    }
00537 
00538    std::vector<double> soln(maxNumEqnsPerNode);
00539 
00540    int numNodes;
00541    CHK_ERR( fei.getNumLocalNodes(numNodes) );
00542 
00543    std::vector<GlobalID> nodes(numNodes);
00544    int* nodesPtr = &nodes[0];
00545 
00546    int checkNumNodes;
00547    CHK_ERR( fei.getLocalNodeIDList( checkNumNodes, nodesPtr, numNodes) );
00548 
00549    if (checkNumNodes != numNodes) {
00550      ERReturn(-1);
00551    }
00552 
00553    FEI_OSTRINGSTREAM fileName;
00554    fileName << solnFileName<<".node."<<solveCounter<<"."<<numProcs<<"."<<localProc;
00555    FEI_OFSTREAM outfile(fileName.str().c_str());
00556 
00557    if (!outfile || outfile.bad()) {
00558       FEI_COUT << "ERROR opening solution output file " << fileName.str() << FEI_ENDL;
00559       return(1);
00560    }
00561 
00562    outfile.setf(IOS_SCIENTIFIC, IOS_FLOATFIELD);
00563 
00564    std::vector<int> offsets(2);
00565 
00566    for(i=0; i<numNodes; ++i) {
00567      CHK_ERR( fei.getNodalSolution(1, &(nodesPtr[i]),
00568            &offsets[0], &soln[0]) );
00569 
00570      int numDOF = offsets[1];
00571 
00572      outfile << nodesPtr[i] << " " << numDOF << FEI_ENDL;
00573      for(int j=0; j<numDOF; j++) {
00574        outfile << soln[j] << " ";
00575      }
00576      outfile << FEI_ENDL;
00577    }
00578 
00579    return(0);
00580 }
00581 
00582 //----------------------------------------------------------------------------
00583 int FEI_tester::save_block_elem_soln(DataReader& data, FEI& fei,
00584              const char* solnFileName,
00585              int numProcs, int localProc,
00586              int solveCounter)
00587 {
00588    int returnValue = 0;
00589    FEI_OSTRINGSTREAM fileName;
00590    fileName << solnFileName<<".elem."<<solveCounter<<"."<<numProcs<<"."<<localProc;
00591    FEI_OFSTREAM outfile(fileName.str().c_str());
00592 
00593    if (!outfile || outfile.bad()) {
00594       FEI_COUT << "ERROR opening elem-solution output file " << fileName.str() << FEI_ENDL;
00595       return(1);
00596    }
00597 
00598    for(int i=0; i<data.numElemBlocks_; i++) {
00599       if (returnValue != 0) break;
00600 
00601       ElemBlock& eb = data.elemBlocks_[i];
00602 
00603       GlobalID blockID = eb.blockID_;
00604       int numElems;
00605       CHK_ERR( fei.getNumBlockElements(blockID, numElems))
00606 
00607       int dofPerElem;
00608       CHK_ERR( fei.getNumBlockElemDOF(blockID, dofPerElem))
00609       int totalNumElemDOF = numElems*dofPerElem;
00610 
00611       if (totalNumElemDOF < 1) {
00612   continue;
00613       }
00614 
00615       GlobalID* elemIDs = new GlobalID[numElems];
00616       if (elemIDs==NULL) return(-1);
00617 
00618       int err = fei.getBlockElemIDList(blockID, numElems, elemIDs);
00619       if (err) returnValue = 1;
00620 
00621       int* offsets = new int[numElems+1];
00622       if (offsets == NULL) return(-1);
00623 
00624       if (totalNumElemDOF > 0) {
00625          double* solnValues = new double[totalNumElemDOF];
00626          if (solnValues == NULL) return(-1);
00627 
00628          err = fei.getBlockElemSolution(blockID, numElems, elemIDs,
00629                                          dofPerElem, solnValues);
00630          if (err) returnValue = 1;
00631 
00632          if (!err) {
00633             for(int j=0; j<numElems; j++) {
00634 
00635                outfile << (int)elemIDs[j] << " " << dofPerElem << FEI_ENDL << "  ";
00636                for(int k=0; k<dofPerElem; k++) {
00637                   outfile << solnValues[j*dofPerElem + k] << " ";
00638                }
00639                outfile << FEI_ENDL;
00640             }
00641          }
00642 
00643          delete [] solnValues;
00644       }
00645 
00646       delete [] elemIDs;
00647       delete [] offsets;
00648    }
00649 
00650    return(0);
00651 }
00652 
00653 //----------------------------------------------------------------------------
00654 int FEI_tester::save_multiplier_soln(DataReader& data, FEI& fei,
00655        const char* solnFileName,
00656                           int numProcs, int localProc, int solveCounter)
00657 {
00658    int numCRs = 0;
00659 
00660    CHK_ERR( fei.getNumCRMultipliers(numCRs) );
00661 
00662    int* globalNumCRs = new int[numProcs];
00663 #ifndef FEI_SER
00664    if (MPI_Allgather(&numCRs, 1, MPI_INT, globalNumCRs, 1, MPI_INT,
00665      comm_) != MPI_SUCCESS) {
00666      ERReturn(-1);
00667    }
00668 #endif
00669 
00670    int localCRStart = 0;
00671    for(int p=0; p<localProc; p++) localCRStart += globalNumCRs[p];
00672 
00673    delete [] globalNumCRs;
00674 
00675    FEI_OSTRINGSTREAM fileName;
00676    fileName << solnFileName<<".mult."<<solveCounter<<"."<<numProcs<<"."<<localProc;
00677    FEI_OFSTREAM outfile(fileName.str().c_str());
00678 
00679    if (!outfile || outfile.bad()) {
00680       FEI_COUT << "ERROR opening mult-solution output file " << fileName.str() << FEI_ENDL;
00681       return(-1);
00682    }
00683 
00684    int* CRIDs = numCRs > 0 ? new int[numCRs] : NULL;
00685    double* results = numCRs > 0 ? new double[numCRs] : NULL;
00686 
00687    if (numCRs > 0 && (CRIDs==NULL || results==NULL)) {
00688      ERReturn(-1);
00689    }
00690 
00691    if (numCRs < 1) {
00692      return(0);
00693    }
00694 
00695    CHK_ERR( fei.getCRMultIDList(numCRs, CRIDs) );
00696 
00697    std::string sname(data_->solverLibraryName_);
00698    if (sname == "FETI") {
00699      for(int ii=0; ii<numCRs; ++ii) results[ii] = -999.99;
00700    }
00701    else {
00702      CHK_ERR( fei.getCRMultipliers(numCRs, CRIDs, results));
00703    }
00704 
00705    for(int i=0; i<numCRs; i++) {
00706       outfile << localCRStart++ << " " << 1 << FEI_ENDL;
00707 
00708       outfile << "   " << results[i] << FEI_ENDL;
00709    }
00710 
00711    delete [] CRIDs;
00712    delete [] results;
00713 
00714    return(0);
00715 }
00716 
00717 //----------------------------------------------------------------------------
00718 int FEI_tester::lsc_matrix_check()
00719 {
00720    if (localProc_ == 0) {
00721      char* current_dir = NULL;
00722      CHK_ERR( fei_test_utils::dirname(data_->solnFileName_.c_str(), current_dir));
00723 
00724      FEI_OSTRINGSTREAM solnMtxName;
00725      solnMtxName<< current_dir<<"/A_TLSC.mtx";
00726      fei::FillableMat solnMtx, checkMtx;
00727      CHK_ERR( SolnCheck::readMatrix(solnMtxName.str().c_str(), numProcs_, solnMtx) );
00728      CHK_ERR( SolnCheck::readMatrix(data_->checkFileName_.c_str(), numProcs_, checkMtx) );
00729      int err = SolnCheck::compareMatrices(solnMtx, checkMtx);
00730      delete [] current_dir;
00731      if (err == 0) {
00732        FEI_COUT << "Utst_fei_lsc: TEST PASSED" << FEI_ENDL;
00733      }
00734      else {
00735        FEI_COUT << "Utst_fei_lsc: TEST FAILED" << FEI_ENDL;
00736        return(-1);
00737      }
00738    }
00739 
00740    return(0);
00741 }

Generated on Wed May 12 21:30:41 2010 for FEI by  doxygen 1.4.7