FEI Version of the Day
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            const int* dof_ids)
00106      {
00107        dbgOut() << "setConnectivity" << FEI_ENDL
00108    << "   elemBlockID: " << elemBlockID << ", elemID: " << elemID 
00109     << ", numNodes: " << numNodes << FEI_ENDL
00110    << "   nodeNumbers: ";
00111        for(int i=0; i<numNodes; i++) {
00112    dbgOut() << nodeNumbers[i] << " ";
00113        }
00114        dbgOut() << FEI_ENDL;
00115 
00116        dbgOut() << "   numDOFPerNode: ";
00117        for(int j=0; j<numNodes; j++) {
00118    dbgOut() << numDofPerNode[j] << " ";
00119        }
00120        dbgOut() << FEI_ENDL;
00121 
00122        return(0);
00123      }
00124 
00125 
00126    int setElemMatrix(int elemBlockID,
00127          int elemID,
00128          int numNodes,
00129          const int* nodeNumbers,
00130          const int* dofPerNode,
00131          const int* dof_ids,
00132          const double *const * coefs)
00133      {
00134        dbgOut() << "setElemMatrix" << FEI_ENDL
00135    << "   elemBlockID: " << elemBlockID << ", elemID: " << elemID << FEI_ENDL
00136    << "   numNodes: " << numNodes << FEI_ENDL;
00137        int i;
00138        dbgOut() << "   nodeNumbers: ";
00139        for(i=0; i<numNodes; i++) {
00140    dbgOut() << nodeNumbers[i] << " ";
00141        }
00142        dbgOut() << FEI_ENDL << "   dofPerNode: ";
00143        int numRows = 0;
00144        for(i=0; i<numNodes; i++) {
00145    dbgOut() << dofPerNode[i] << " ";
00146    numRows += dofPerNode[i];
00147        }
00148        dbgOut() << FEI_ENDL << "   coefs:" << FEI_ENDL;
00149        for(i=0; i<numRows; i++) {
00150    dbgOut() << "      ";
00151    for(int j=0; j<numRows; j++) {
00152      dbgOut() << coefs[i][j] << " ";
00153    }
00154    dbgOut() << FEI_ENDL;
00155        }
00156 
00157        return(0);
00158      }
00159 
00160 
00161    int setElemVector(int elemBlockID,
00162          int elemID,
00163          int numNodes,
00164          const int* nodeNumbers,
00165          const int* dofPerNode,
00166          const int* dof_ids,
00167          const double* coefs)
00168      {
00169        dbgOut() << "setElemVector" << FEI_ENDL
00170    << "   elemBlockID: " << elemBlockID << ", elemID: " << elemID << FEI_ENDL
00171    << "   numNodes: " << numNodes << FEI_ENDL;
00172        int i;
00173        dbgOut() << "   nodeNumbers: ";
00174        for(i=0; i<numNodes; i++) {
00175    dbgOut() << nodeNumbers[i] << " ";
00176        }
00177        dbgOut() << FEI_ENDL << "   dofPerNode: ";
00178        int numRows = 0;
00179        for(i=0; i<numNodes; i++) {
00180    dbgOut() << dofPerNode[i] << " ";
00181    numRows += dofPerNode[i];
00182        }
00183        dbgOut() << FEI_ENDL << "   coefs:" << FEI_ENDL << "      ";
00184        for(i=0; i<numRows; i++) {
00185    dbgOut() << coefs[i] << " ";
00186        }
00187        dbgOut() << FEI_ENDL;
00188        return(0);
00189      }
00190 
00191    int setDirichletBCs(int numBCs,
00192            const int* nodeNumbers,
00193            const int* dofOffsets,
00194            const double* values)
00195      {
00196        dbgOut() << "setDirichletBCs" << FEI_ENDL
00197    << "   numBCs: " << numBCs << FEI_ENDL;
00198        for(int i=0; i<numBCs; i++) {
00199    dbgOut() << "     nodeNumber: " << nodeNumbers[i] << ", "
00200      << "dof-offset: " << dofOffsets[i] << ", value: " << values[i]<<FEI_ENDL;
00201        }
00202 
00203        return(0);
00204      }
00205 
00206    int sumIntoMatrix(int numRowNodes,
00207          const int* rowNodeNumbers,
00208          const int* rowDofOffsets,
00209          const int* numColNodesPerRow,
00210          const int* colNodeNumbers,
00211          const int* colDofOffsets,
00212          const double* coefs)
00213      {
00214        dbgOut() << "sumIntoMatrix, numRowNodes: " << numRowNodes << FEI_ENDL;
00215        int offset = 0;
00216        for(int i=0; i<numRowNodes; i++) {
00217    dbgOut() << "   rowNodeNumber " << rowNodeNumbers[i]
00218      << ", rowDofOffset " << rowDofOffsets[i] << FEI_ENDL;
00219    for(int j=0; j<numColNodesPerRow[i]; j++) {
00220      dbgOut() << "      colNodeNumber " << colNodeNumbers[offset]
00221         << ", colDofOffset " << colDofOffsets[offset]
00222         << ", value: " << coefs[offset]<<FEI_ENDL;
00223      offset++;
00224    }
00225        }
00226 
00227        return(0);
00228      }
00229 
00230    int sumIntoRHSVector(int numNodes,
00231          const int* nodeNumbers,
00232          const int* dofOffsets,
00233          const double* coefs)
00234      {
00235        dbgOut() << "sumIntoRHSVector, numNodes: " << numNodes << FEI_ENDL;
00236        for(int i=0; i<numNodes; i++) {
00237    dbgOut() << "   nodeNumber " << nodeNumbers[i]
00238      << ", dof-offset " << dofOffsets[i] << ", value: " << coefs[i]<<FEI_ENDL;
00239        }
00240 
00241        return(0);
00242      }
00243 
00244    int putIntoRHSVector(int numNodes,
00245          const int* nodeNumbers,
00246          const int* dofOffsets,
00247          const double* coefs)
00248      {
00249        dbgOut() << "putIntoRHSVector, numNodes: " << numNodes << FEI_ENDL;
00250        for(int i=0; i<numNodes; i++) {
00251    dbgOut() << "   nodeNumber " << nodeNumbers[i]
00252      << ", dof-offset " << dofOffsets[i] << ", value: " << coefs[i]<<FEI_ENDL;
00253        }
00254 
00255        return(0);
00256      }
00257 
00258    int loadComplete()
00259      {
00260        dbgOut() << "loadComplete" << FEI_ENDL;
00261 
00262        return(0);
00263      }
00264 
00274    int launchSolver(int& solveStatus, int& iterations)
00275      {
00276        dbgOut() << "launchSolver" << FEI_ENDL;
00277 
00278        solveStatus = 0;
00279        iterations = 0;
00280 
00281        return(0);
00282      }
00283 
00284    int reset()
00285      {
00286        dbgOut() << "reset" << FEI_ENDL;
00287        return(0);
00288      }
00289 
00290    int resetRHSVector()
00291      {
00292        dbgOut() << "resetRHSVector" << FEI_ENDL;
00293        return(0);
00294      }
00295 
00296    int resetMatrix()
00297      {
00298        dbgOut() << "resetMatrix" << FEI_ENDL;
00299        return(0);
00300      }
00301 
00302    int deleteConstraints()
00303      {
00304        dbgOut() << "deleteConstraints" << FEI_ENDL;
00305        return(0);
00306      }
00307 
00308    int getSolnEntry(int nodeNumber,
00309         int dofOffset,
00310         double& value)
00311      {
00312        dbgOut() << "getSolnEntry, nodeNumber: " << nodeNumber 
00313    << ", dofOffset: " << dofOffset << FEI_ENDL;
00314 
00315        value = -999.99;
00316 
00317        return(0);
00318     }
00319 
00320    int getMultiplierSoln(int CRID, double& lagrangeMultiplier)
00321      {
00322        lagrangeMultiplier = -999.99;
00323        return(0);
00324      }
00325 
00337    int putNodalFieldData(int fieldID,
00338        int fieldSize,
00339        int numNodes,
00340        const int* nodeNumbers,
00341        const double* coefs)
00342      {
00343        dbgOut() << "putNodalFieldData, fieldID: " << fieldID << ", fieldSize: "
00344    << fieldSize << FEI_ENDL;
00345        int offset = 0;
00346        for(int i=0; i<numNodes; i++) {
00347    dbgOut() << "   nodeNumber " << nodeNumbers[i] << ", coefs: ";
00348    for(int j=0; j<fieldSize; j++) {
00349      dbgOut() << coefs[offset++] << " ";
00350    }
00351    dbgOut() << FEI_ENDL;
00352        }
00353        
00354        return(0);
00355      }
00356 
00357    int setMultiplierCR(int CRID,
00358            int numNodes,
00359            const int* nodeNumbers,
00360            const int* dofOffsets,
00361            const double* coefWeights,
00362            double rhsValue)
00363      {
00364        dbgOut() << "setMultiplierCR, CRID: " << CRID << ", numNodes: " << numNodes << FEI_ENDL;
00365        for(int i=0; i<numNodes; i++) {
00366    dbgOut() << "   nodeNumber " << nodeNumbers[i] << ", dof-offset "
00367       << dofOffsets[i] << ", coefWeight: " << coefWeights[i] <<FEI_ENDL;
00368        }
00369 
00370        dbgOut() << "   rhsValue: " << rhsValue << FEI_ENDL;
00371 
00372        return(0);
00373      }
00374 
00375    int setPenaltyCR(int CRID,
00376         int numNodes,
00377         const int* nodeNumbers,
00378         const int* dofOffsets,
00379         const double* coefWeights,
00380         double penaltyValue,
00381         double rhsValue)
00382      {
00383        dbgOut() << "setPenaltyCR, CRID: " << CRID << ", numNodes: " << numNodes << FEI_ENDL;
00384        for(int i=0; i<numNodes; i++) {
00385    dbgOut() << "   nodeNumber " << nodeNumbers[i] << ", dof-offset "
00386       << dofOffsets[i] << ", coefWeight: " << coefWeights[i] <<FEI_ENDL;
00387        }
00388 
00389        dbgOut() << "   penaltyValue: " << penaltyValue << FEI_ENDL;
00390        dbgOut() << "   rhsValue: " << rhsValue << FEI_ENDL;
00391 
00392        return(0);
00393      }
00394 
00395  private:
00396    int setDebugLog(int debugOutputLevel, const char* path);
00397 
00398    MPI_Comm comm_;
00399    int numProcs_, localProc_;
00400 
00401    int debugOutputLevel_;
00402    char* dbgPath_;
00403    FEI_OSTREAM* dbgOStreamPtr_;
00404    bool dbgFileOpened_;
00405    FEI_OFSTREAM* dbgFStreamPtr_;
00406 };
00407 
00408 #endif // _FEData_h_
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends