FEData.hpp

00001 #ifndef _FEData_h_
00002 #define _FEData_h_
00003 
00004 /*--------------------------------------------------------------------*/
00005 /*    Copyright 2005 Sandia Corporation.                              */
00006 /*    Under the terms of Contract DE-AC04-94AL85000, there is a       */
00007 /*    non-exclusive license for use of this work by or on behalf      */
00008 /*    of the U.S. Government.  Export of this program may require     */
00009 /*    a license from the United States Government.                    */
00010 /*--------------------------------------------------------------------*/
00011 
00012 #include <fei_iostream.hpp>
00013 #include <fei_fstream.hpp>
00014 
00015 #include <fei_mpi.h>
00016 
00017 #define dbgOut() if (debugOutputLevel_ > 0) *dbgOStreamPtr_
00018 
00019 #include <fei_defs.h>
00020 #include <fei_FiniteElementData.hpp>
00021 
00022 #include <snl_fei_Utils.hpp>
00023 
00028 class FEData : public virtual FiniteElementData {
00029  public:
00031   FEData(MPI_Comm comm)
00032     :
00033     comm_(comm), numProcs_(1), localProc_(0),
00034     debugOutputLevel_(0),
00035     dbgPath_(NULL),
00036     dbgOStreamPtr_(NULL),
00037     dbgFileOpened_(false)
00038     {
00039 #ifndef FEI_SER
00040       if (MPI_Comm_rank(comm_, &localProc_) != MPI_SUCCESS) MPI_Abort(comm_,-1);
00041       if (MPI_Comm_size(comm_, &numProcs_) != MPI_SUCCESS) MPI_Abort(comm_,-1);
00042 #endif
00043       setDebugLog(0, ".");
00044     }
00045 
00046   virtual ~FEData()
00047     {
00048       if (dbgFileOpened_ == true) { dbgFStreamPtr_->close(); }
00049 
00050       delete [] dbgPath_;
00051       delete dbgOStreamPtr_;
00052     }
00053 
00059   int parameters(int numParams, char** params);
00060 
00068   int setLookup(Lookup& lookup)
00069     {
00070       dbgOut() << "setLookup" << FEI_ENDL;
00071       return(0);
00072     }
00073 
00074 
00078   int describeStructure(int numElemBlocks,
00079       const int* numElemsPerBlock,
00080       const int* numNodesPerElem,
00081       const int* elemMatrixSizePerBlock,
00082       int totalNumNodes,
00083       int numSharedNodes,
00084       int numMultCRs)
00085     {
00086       dbgOut() << "describeStructure" << FEI_ENDL
00087   << "   numElemBlocks: " << numElemBlocks << FEI_ENDL;
00088       for(int i=0; i<numElemBlocks; i++) {
00089   dbgOut() << "   elem-block " << i << ": " << FEI_ENDL
00090     << "      number of elements: " << numElemsPerBlock[i] << FEI_ENDL
00091     << "      nodes per element:  " << numNodesPerElem[i] << FEI_ENDL;
00092   dbgOut() << "         elemMatrixSizePerBlock: "
00093       << elemMatrixSizePerBlock[i] << FEI_ENDL;
00094       }
00095       return(0);
00096     }
00097 
00100    int setConnectivity(int elemBlockID,
00101            int elemID,
00102            int numNodes,
00103            const int* nodeNumbers,
00104            const int* numDofPerNode)
00105      {
00106        dbgOut() << "setConnectivity" << FEI_ENDL
00107    << "   elemBlockID: " << elemBlockID << ", elemID: " << elemID 
00108     << ", numNodes: " << numNodes << FEI_ENDL
00109    << "   nodeNumbers: ";
00110        for(int i=0; i<numNodes; i++) {
00111    dbgOut() << nodeNumbers[i] << " ";
00112        }
00113        dbgOut() << FEI_ENDL;
00114 
00115        dbgOut() << "   numDOFPerNode: ";
00116        for(int j=0; j<numNodes; j++) {
00117    dbgOut() << numDofPerNode[j] << " ";
00118        }
00119        dbgOut() << FEI_ENDL;
00120 
00121        return(0);
00122      }
00123 
00124 
00125    int setElemMatrix(int elemBlockID,
00126          int elemID,
00127          int numNodes,
00128          const int* nodeNumbers,
00129          const int* dofPerNode,
00130          const double *const * coefs)
00131      {
00132        dbgOut() << "setElemMatrix" << FEI_ENDL
00133    << "   elemBlockID: " << elemBlockID << ", elemID: " << elemID << FEI_ENDL
00134    << "   numNodes: " << numNodes << FEI_ENDL;
00135        int i;
00136        dbgOut() << "   nodeNumbers: ";
00137        for(i=0; i<numNodes; i++) {
00138    dbgOut() << nodeNumbers[i] << " ";
00139        }
00140        dbgOut() << FEI_ENDL << "   dofPerNode: ";
00141        int numRows = 0;
00142        for(i=0; i<numNodes; i++) {
00143    dbgOut() << dofPerNode[i] << " ";
00144    numRows += dofPerNode[i];
00145        }
00146        dbgOut() << FEI_ENDL << "   coefs:" << FEI_ENDL;
00147        for(i=0; i<numRows; i++) {
00148    dbgOut() << "      ";
00149    for(int j=0; j<numRows; j++) {
00150      dbgOut() << coefs[i][j] << " ";
00151    }
00152    dbgOut() << FEI_ENDL;
00153        }
00154 
00155        return(0);
00156      }
00157 
00158 
00159    int setElemVector(int elemBlockID,
00160          int elemID,
00161          int numNodes,
00162          const int* nodeNumbers,
00163          const int* dofPerNode,
00164          const double* coefs)
00165      {
00166        dbgOut() << "setElemVector" << FEI_ENDL
00167    << "   elemBlockID: " << elemBlockID << ", elemID: " << elemID << FEI_ENDL
00168    << "   numNodes: " << numNodes << FEI_ENDL;
00169        int i;
00170        dbgOut() << "   nodeNumbers: ";
00171        for(i=0; i<numNodes; i++) {
00172    dbgOut() << nodeNumbers[i] << " ";
00173        }
00174        dbgOut() << FEI_ENDL << "   dofPerNode: ";
00175        int numRows = 0;
00176        for(i=0; i<numNodes; i++) {
00177    dbgOut() << dofPerNode[i] << " ";
00178    numRows += dofPerNode[i];
00179        }
00180        dbgOut() << FEI_ENDL << "   coefs:" << FEI_ENDL << "      ";
00181        for(i=0; i<numRows; i++) {
00182    dbgOut() << coefs[i] << " ";
00183        }
00184        dbgOut() << FEI_ENDL;
00185        return(0);
00186      }
00187 
00188    int setDirichletBCs(int numBCs,
00189            const int* nodeNumbers,
00190            const int* dofOffsets,
00191            const double* values)
00192      {
00193        dbgOut() << "setDirichletBCs" << FEI_ENDL
00194    << "   numBCs: " << numBCs << FEI_ENDL;
00195        for(int i=0; i<numBCs; i++) {
00196    dbgOut() << "     nodeNumber: " << nodeNumbers[i] << ", "
00197      << "dof-offset: " << dofOffsets[i] << ", value: " << values[i]<<FEI_ENDL;
00198        }
00199 
00200        return(0);
00201      }
00202 
00203    int sumIntoMatrix(int numRowNodes,
00204          const int* rowNodeNumbers,
00205          const int* rowDofOffsets,
00206          const int* numColNodesPerRow,
00207          const int* colNodeNumbers,
00208          const int* colDofOffsets,
00209          const double* coefs)
00210      {
00211        dbgOut() << "sumIntoMatrix, numRowNodes: " << numRowNodes << FEI_ENDL;
00212        int offset = 0;
00213        for(int i=0; i<numRowNodes; i++) {
00214    dbgOut() << "   rowNodeNumber " << rowNodeNumbers[i]
00215      << ", rowDofOffset " << rowDofOffsets[i] << FEI_ENDL;
00216    for(int j=0; j<numColNodesPerRow[i]; j++) {
00217      dbgOut() << "      colNodeNumber " << colNodeNumbers[offset]
00218         << ", colDofOffset " << colDofOffsets[offset]
00219         << ", value: " << coefs[offset]<<FEI_ENDL;
00220      offset++;
00221    }
00222        }
00223 
00224        return(0);
00225      }
00226 
00227    int sumIntoRHSVector(int numNodes,
00228          const int* nodeNumbers,
00229          const int* dofOffsets,
00230          const double* coefs)
00231      {
00232        dbgOut() << "sumIntoRHSVector, numNodes: " << numNodes << FEI_ENDL;
00233        for(int i=0; i<numNodes; i++) {
00234    dbgOut() << "   nodeNumber " << nodeNumbers[i]
00235      << ", dof-offset " << dofOffsets[i] << ", value: " << coefs[i]<<FEI_ENDL;
00236        }
00237 
00238        return(0);
00239      }
00240 
00241    int putIntoRHSVector(int numNodes,
00242          const int* nodeNumbers,
00243          const int* dofOffsets,
00244          const double* coefs)
00245      {
00246        dbgOut() << "putIntoRHSVector, numNodes: " << numNodes << FEI_ENDL;
00247        for(int i=0; i<numNodes; i++) {
00248    dbgOut() << "   nodeNumber " << nodeNumbers[i]
00249      << ", dof-offset " << dofOffsets[i] << ", value: " << coefs[i]<<FEI_ENDL;
00250        }
00251 
00252        return(0);
00253      }
00254 
00255    int loadComplete()
00256      {
00257        dbgOut() << "loadComplete" << FEI_ENDL;
00258 
00259        return(0);
00260      }
00261 
00271    int launchSolver(int& solveStatus, int& iterations)
00272      {
00273        dbgOut() << "launchSolver" << FEI_ENDL;
00274 
00275        solveStatus = 0;
00276        iterations = 0;
00277 
00278        return(0);
00279      }
00280 
00281    int reset()
00282      {
00283        dbgOut() << "reset" << FEI_ENDL;
00284        return(0);
00285      }
00286 
00287    int resetRHSVector()
00288      {
00289        dbgOut() << "resetRHSVector" << FEI_ENDL;
00290        return(0);
00291      }
00292 
00293    int resetMatrix()
00294      {
00295        dbgOut() << "resetMatrix" << FEI_ENDL;
00296        return(0);
00297      }
00298 
00299    int deleteConstraints()
00300      {
00301        dbgOut() << "deleteConstraints" << FEI_ENDL;
00302        return(0);
00303      }
00304 
00305    int getSolnEntry(int nodeNumber,
00306         int dofOffset,
00307         double& value)
00308      {
00309        dbgOut() << "getSolnEntry, nodeNumber: " << nodeNumber 
00310    << ", dofOffset: " << dofOffset << FEI_ENDL;
00311 
00312        value = -999.99;
00313 
00314        return(0);
00315     }
00316 
00317    int getMultiplierSoln(int CRID, double& lagrangeMultiplier)
00318      {
00319        lagrangeMultiplier = -999.99;
00320        return(0);
00321      }
00322 
00334    int putNodalFieldData(int fieldID,
00335        int fieldSize,
00336        int numNodes,
00337        const int* nodeNumbers,
00338        const double* coefs)
00339      {
00340        dbgOut() << "putNodalFieldData, fieldID: " << fieldID << ", fieldSize: "
00341    << fieldSize << FEI_ENDL;
00342        int offset = 0;
00343        for(int i=0; i<numNodes; i++) {
00344    dbgOut() << "   nodeNumber " << nodeNumbers[i] << ", coefs: ";
00345    for(int j=0; j<fieldSize; j++) {
00346      dbgOut() << coefs[offset++] << " ";
00347    }
00348    dbgOut() << FEI_ENDL;
00349        }
00350        
00351        return(0);
00352      }
00353 
00354    int setMultiplierCR(int CRID,
00355            int numNodes,
00356            const int* nodeNumbers,
00357            const int* dofOffsets,
00358            const double* coefWeights,
00359            double rhsValue)
00360      {
00361        dbgOut() << "setMultiplierCR, CRID: " << CRID << ", numNodes: " << numNodes << FEI_ENDL;
00362        for(int i=0; i<numNodes; i++) {
00363    dbgOut() << "   nodeNumber " << nodeNumbers[i] << ", dof-offset "
00364       << dofOffsets[i] << ", coefWeight: " << coefWeights[i] <<FEI_ENDL;
00365        }
00366 
00367        dbgOut() << "   rhsValue: " << rhsValue << FEI_ENDL;
00368 
00369        return(0);
00370      }
00371 
00372    int setPenaltyCR(int CRID,
00373         int numNodes,
00374         const int* nodeNumbers,
00375         const int* dofOffsets,
00376         const double* coefWeights,
00377         double penaltyValue,
00378         double rhsValue)
00379      {
00380        dbgOut() << "setPenaltyCR, CRID: " << CRID << ", numNodes: " << numNodes << FEI_ENDL;
00381        for(int i=0; i<numNodes; i++) {
00382    dbgOut() << "   nodeNumber " << nodeNumbers[i] << ", dof-offset "
00383       << dofOffsets[i] << ", coefWeight: " << coefWeights[i] <<FEI_ENDL;
00384        }
00385 
00386        dbgOut() << "   penaltyValue: " << penaltyValue << FEI_ENDL;
00387        dbgOut() << "   rhsValue: " << rhsValue << FEI_ENDL;
00388 
00389        return(0);
00390      }
00391 
00392  private:
00393    int setDebugLog(int debugOutputLevel, const char* path);
00394 
00395    MPI_Comm comm_;
00396    int numProcs_, localProc_;
00397 
00398    int debugOutputLevel_;
00399    char* dbgPath_;
00400    FEI_OSTREAM* dbgOStreamPtr_;
00401    bool dbgFileOpened_;
00402    FEI_OFSTREAM* dbgFStreamPtr_;
00403 };
00404 
00405 #endif // _FEData_h_

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