SNL_FEI_Structure.hpp

00001 #ifndef _SNL_FEI_Structure_hpp_
00002 #define _SNL_FEI_Structure_hpp_
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_fwd.hpp>
00013 #include <fei_defs.h>
00014 #include "fei_TemplateUtils.hpp"
00015 #include <snl_fei_PointBlockMap.hpp>
00016 #include <fei_EqnBuffer.hpp>
00017 #include <fei_CSRMat.hpp>
00018 #include <fei_CSVec.hpp>
00019 
00020 #include <fei_NodeCommMgr.hpp>
00021 #include <fei_NodeDatabase.hpp>
00022 #include <fei_NodeDescriptor.hpp>
00023 
00024 #include <fei_Lookup.hpp>
00025 
00041 class SNL_FEI_Structure : public Lookup {
00042  public:
00051    SNL_FEI_Structure(MPI_Comm comm);
00052 
00054    virtual ~SNL_FEI_Structure();
00055 
00071    int parameters(int numParams, const char*const* paramStrings);
00072 
00073    int initFields(int numFields, const int* fieldSizes, const int* fieldIDs);
00074 
00075    int initElemBlock(GlobalID elemBlockID,
00076          int numElements,
00077          int numNodesPerElement,
00078          const int* numFieldsPerNode,
00079          const int* const* nodalFieldIDs,
00080          int numElemDofFieldsPerElement,
00081          const int* elemDofFieldIDs,
00082          int interleaveStrategy);
00083 
00084    int initElem(GlobalID elemBlockID,
00085     GlobalID elemID,
00086     const GlobalID* elemConn);
00087 
00088    int initSlaveVariable(GlobalID slaveNodeID, 
00089        int slaveFieldID,
00090        int offsetIntoSlaveField,
00091        int numMasterNodes,
00092        const GlobalID* masterNodeIDs,
00093        const int* masterFieldIDs,
00094        const double* weights,
00095        double rhsValue);
00096 
00097    int deleteMultCRs();
00098 
00099    int initSharedNodes(int numSharedNodes,
00100            const GlobalID *sharedNodeIDs,  
00101            const int* numProcsPerNode, 
00102            const int *const *sharingProcIDs);
00103 
00104    // constraint relation initialization
00105    //- lagrange multiplier formulation
00106    int initCRMult(int numCRNodes,
00107       const GlobalID* CRNodes,
00108       const int *CRFields,
00109       int& CRID);
00110 
00111    // - penalty function formulation
00112    int initCRPen(int numCRNodes,
00113      const GlobalID* CRNodes, 
00114      const int *CRFields,
00115      int& CRID); 
00116 
00117    int initComplete(bool generateGraph = true);
00118 
00119    std::map<int,int>& getFieldDatabase() { return(*fieldDatabase_); };
00120 
00122    const int* getFieldIDsPtr()
00123      {
00124        int len = fieldDatabase_->size();
00125        workarray_.resize(len*2);
00126        fei::copyToArrays<std::map<int,int> >(*fieldDatabase_, len,
00127              &workarray_[0],
00128              &workarray_[0]+len);
00129        return( &workarray_[0] );
00130      }
00131 
00133    const int* getFieldSizesPtr()
00134      {
00135        int len = fieldDatabase_->size();
00136        workarray_.resize(len*2);
00137        fei::copyToArrays<std::map<int,int> >(*fieldDatabase_, len,
00138              &workarray_[0],
00139              &workarray_[0]+len);
00140        return( &workarray_[0]+len );
00141      }
00142 
00144    int getNumFields() { return( fieldDatabase_->size() ); };
00145 
00147    int getFieldSize(int fieldID)
00148      {
00149        std::map<int,int>::const_iterator
00150    f_iter = fieldDatabase_->find(fieldID);
00151 
00152        return(f_iter != fieldDatabase_->end() ? (*f_iter).second : -1);
00153      }
00154 
00156    bool isInLocalElement(int nodeNumber);
00157 
00158    const int* getNumFieldsPerNode(GlobalID blockID);
00159 
00160    const int* const* getFieldIDsTable(GlobalID blockID);
00161 
00166    int getEqnNumber(int nodeNumber, int fieldID);
00167 
00174    int getOwnerProcForEqn(int eqn);
00175 
00176 
00178    //now the element-block functions
00179 
00180    int getNumElemBlocks() {return(blockIDs_.size());};
00181    const GlobalID* getElemBlockIDs() {return(&blockIDs_[0]);};
00182 
00183    void getElemBlockInfo(GlobalID blockID,
00184                          int& interleaveStrategy, int& lumpingStrategy,
00185                          int& numElemDOF, int& numElements,
00186                          int& numNodesPerElem, int& numEqnsPerElem);
00187 
00188    int addBlock(GlobalID blockID);
00189 
00190    int getBlockDescriptor(GlobalID blockID, BlockDescriptor*& block);
00191 
00195    int getBlockDescriptor_index(int index, BlockDescriptor*& block);
00196 
00198    int getIndexOfBlock(GlobalID blockID);
00199 
00200    int allocateBlockConnectivity(GlobalID blockID);
00201    void destroyConnectivityTables();
00202 
00203    ConnectivityTable& getBlockConnectivity(GlobalID blockID);
00204 
00205    void getScatterIndices_ID(GlobalID blockID, GlobalID elemID,
00206                              int interleaveStrategy,
00207                              int* scatterIndices);
00208 
00209    void getScatterIndices_index(int blockIndex, int elemIndex,
00210                              int interleaveStrategy,
00211                              int* scatterIndices);
00212 
00213    int getBlkScatterIndices_index(int blockIndex,
00214           int elemIndex,
00215           int* scatterIndices);
00216 
00217    void getScatterIndices_ID(GlobalID blockID, GlobalID elemID,
00218                              int interleaveStrategy,
00219                              int* scatterIndices,
00220            int* blkScatterIndices,
00221            int* blkSizes);
00222 
00223    void getScatterIndices_index(int blockIndex, int elemIndex,
00224                              int interleaveStrategy,
00225                              int* scatterIndices,
00226            int* blkScatterIndices,
00227         int* blkSizes);
00228 
00229 
00231    //now the shared-node lookup functions from the Lookup interface.
00232 
00233    int getNumSharedNodes() {return(nodeCommMgr_->getNumSharedNodes());};
00234 
00235    const int* getSharedNodeNumbers() {
00236       return(&(nodeCommMgr_->getSharedNodeNumbers())[0]);
00237    };
00238 
00239    const int* getSharedNodeProcs(int nodeNumber) {
00240       int index = nodeCommMgr_->getSharedNodeIndex_num(nodeNumber);
00241       if (index < 0) return(NULL);
00242       return(&(nodeCommMgr_->getSharedNodeProcs(index))[0]);
00243    };
00244 
00245    int getNumSharingProcs(int nodeNumber) {
00246       int index = nodeCommMgr_->getSharedNodeIndex_num(nodeNumber);
00247       if (index < 0) return(-1);
00248       return(nodeCommMgr_->getSharedNodeProcs(index).size());
00249    };
00250 
00251    int getNumSubdomains(int nodeNumber) {
00252      const NodeDescriptor* node = NULL;
00253      int err = nodeDatabase_->getNodeWithNumber(nodeNumber, node);
00254      if (err != 0) return(-1);
00255      GlobalID nodeID = node->getGlobalNodeID();
00256      return(nodeCommMgr_->getSharedNodeNumSubdomains(nodeID));
00257    }
00258 
00259    int* getSubdomainList(int nodeNumber) {
00260      const NodeDescriptor* node = NULL;
00261      int err = nodeDatabase_->getNodeWithNumber(nodeNumber, node);
00262      if (err != 0) return(NULL);
00263      GlobalID nodeID = node->getGlobalNodeID();
00264      return(&(*(nodeCommMgr_->getSharedNodeSubdomainList(nodeID)))[0]);
00265 
00266    }
00267 
00268    int translateToReducedNodeNumber(int nodeNumber, int proc);
00269 
00271 
00273    int getAssociatedNodeNumber(int eqnNumber)
00274      {
00275        int eqn = translateFromReducedEqn(eqnNumber);
00276        int nodeNumber = nodeDatabase_->getAssociatedNodeNumber(eqn);
00277        int reducedNodeNumber = -1;
00278        if (nodeNumber >= 0) {
00279    reducedNodeNumber = translateToReducedNodeNumber(nodeNumber, localProc_);
00280        }
00281        return( reducedNodeNumber );
00282      }
00283 
00285    int getAssociatedFieldID(int eqnNumber)
00286      {
00287        int eqn = translateFromReducedEqn(eqnNumber);
00288        return( nodeDatabase_->getAssociatedFieldID(eqn) );
00289      }
00290 
00291 
00293    //now the point-eqn to block-eqn queries...
00294 
00295    bool isExactlyBlkEqn(int ptEqn) {
00296       return(blkEqnMapper_->isExactlyBlkEqn(ptEqn));
00297    };
00298 
00299    int ptEqnToBlkEqn(int ptEqn) {
00300       return(blkEqnMapper_->eqnToBlkEqn(ptEqn));
00301    };
00302  
00303    int getOffsetIntoBlkEqn(int blkEqn, int ptEqn) {
00304       return(blkEqnMapper_->getBlkEqnOffset(blkEqn, ptEqn));
00305    };
00306 
00307    int getBlkEqnSize(int blkEqn) {
00308       return(blkEqnMapper_->getBlkEqnSize(blkEqn));
00309    }
00311 
00312 
00313    int getNumActiveNodes() {return(nodeDatabase_->getNodeIDs().size());}
00314 
00315    NodeDatabase& getNodeDatabase() { return( *nodeDatabase_ ); }
00316 
00317    std::map<GlobalID,int>& getActiveNodeIDList()
00318      { return( nodeDatabase_->getNodeIDs() ); }
00319 
00320    std::vector<int>& getGlobalNodeOffsets() {return(globalNodeOffsets_);}
00321    std::vector<int>& getGlobalEqnOffsets() {return(globalEqnOffsets_);}
00322    std::vector<int>& getGlobalBlkEqnOffsets() {return(globalBlkEqnOffsets_);}
00323 
00324    NodeCommMgr& getNodeCommMgr() {return(*nodeCommMgr_);}
00325    EqnCommMgr&  getEqnCommMgr()  {return(*eqnCommMgr_ );}
00326 
00327    void initializeEqnCommMgr();
00328 
00329    void getEqnInfo(int& numGlobalEqns, int& numLocalEqns,
00330                    int& localStartRow, int& localEndRow);
00331 
00332    int getEqnNumbers(GlobalID ID, int idType, int fieldID,
00333          int& numEqns, int* eqnNumbers);
00334 
00335    int getEqnNumbers(int numIDs, const GlobalID* IDs,
00336          int idType, int fieldID,
00337          int& numEqns, int* eqnNumbers);
00338 
00339    void getEqnBlkInfo(int& numGlobalEqnBlks, int& numLocalEqnBlks,
00340                       int& localBlkOffset);
00341 
00342    snl_fei::PointBlockMap& getBlkEqnMapper() {return(*blkEqnMapper_);}
00343 
00344    void destroyMatIndices();
00345 
00346    int getNumMultConstRecords() {return(multCRs_.size());};
00347 
00348    std::map<GlobalID,snl_fei::Constraint<GlobalID>*>&
00349      getMultConstRecords()
00350      {return(multCRs_);};
00351 
00352    int getMultConstRecord(int CRID, snl_fei::Constraint<GlobalID>*& multCR)
00353      {
00354        std::map<int,snl_fei::Constraint<GlobalID>*>::iterator
00355    cr_iter = multCRs_.find(CRID);
00356        int returncode = -1;
00357        if (cr_iter != multCRs_.end()) {
00358    multCR = (*cr_iter).second;
00359    returncode = 0;
00360        }
00361 
00362        return( returncode );
00363      }
00364 
00365    int getNumPenConstRecords() {return(penCRs_.size());}
00366    std::map<GlobalID,snl_fei::Constraint<GlobalID>*>&
00367      getPenConstRecords()
00368      { return(penCRs_); }
00369 
00370    int getPenConstRecord(int CRID, snl_fei::Constraint<GlobalID>*& penCR)
00371      {
00372        std::map<int,snl_fei::Constraint<GlobalID>*>::iterator
00373    cr_iter = penCRs_.find(CRID);
00374        int returncode = -1;
00375        if (cr_iter != penCRs_.end()) {
00376    penCR = (*cr_iter).second;
00377    returncode = 0;
00378        }
00379 
00380        return( returncode );
00381      }
00382 
00383    void addSlaveVariable(SlaveVariable* svar) {slaveVars_->push_back(svar);}
00384 
00385    int calculateSlaveEqns(MPI_Comm comm);
00386 
00387    fei::FillableMat* getSlaveDependencies() {return(slaveMatrix_);}
00388 
00389    EqnBuffer* getSlaveEqns() { return(slaveEqns_); }
00390 
00391    int numSlaveEquations() { return(numSlvs_); }
00392 
00395    bool isSlaveEqn(int eqn);
00396 
00405    bool translateToReducedEqn(int eqn, int& reducedEqn);
00406 
00411    int translateToReducedEqns(EqnCommMgr& eqnCommMgr);
00412 
00417    int translateToReducedEqns(EqnBuffer& eqnBuf);
00418 
00423    int translateToReducedEqns(ProcEqns& procEqns);
00424 
00429    int translateMatToReducedEqns(fei::CSRMat& mat);
00430 
00436    int translateFromReducedEqn(int reducedEqn);
00437 
00444    int getMasterEqnNumbers(int slaveEqn, std::vector<int>*& masterEqns);
00445 
00452    int getMasterEqnCoefs(int slaveEqn, std::vector<double>*& masterCoefs);
00453 
00461    int getMasterEqnRHS(int slaveEqn, double& rhsValue);
00462 
00463    int getNumGlobalEqns() { return( numGlobalEqns_ ); }
00464    int getNumLocalEqns() { return( numLocalEqns_ ); }
00465    int getFirstLocalEqn() { return( localStartRow_ ); }
00466    int getLastLocalEqn() { return( localEndRow_ ); }
00467 
00468    int getFirstReducedEqn() { return( reducedStartRow_ ); }
00469    int getLastReducedEqn()  { return( reducedEndRow_   ); }
00470 
00471    int getNumGlobalEqnBlks() { return( numGlobalEqnBlks_ ); }
00472    int getNumLocalEqnBlks() { return( numLocalEqnBlks_ ); }
00473    int getNumLocalReducedEqnBlks() { return( numLocalReducedEqnBlks_ ); }
00474    int getGlobalMaxBlkSize() { return(globalMaxBlkSize_); }
00475 
00476    int getNumLocalReducedEqns() { return( numLocalReducedRows_ ); }
00477 
00478    int getMatrixRowLengths(std::vector<int>& rowLengths);
00479    int getMatrixStructure(int** colIndices, std::vector<int>& rowLengths);
00480 
00481    int getMatrixStructure(int** ptColIndices, std::vector<int>& ptRowLengths,
00482         int** blkColIndices, int* blkIndices_1D,
00483         std::vector<int>& blkRowLengths,
00484         std::vector<int>& numPtRowsPerBlkRow);
00485 
00486    static int gatherSlaveEqns(MPI_Comm comm,
00487             EqnCommMgr* eqnCommMgr,
00488             EqnBuffer* slaveEqns);
00489 
00490    static int removeCouplings(EqnBuffer& eqnbuf, int& levelsOfCoupling);
00491 
00492    int calcTotalNumElemDOF();
00493    int calcNumMultCREqns();
00494 
00495    MPI_Comm getCommunicator() const { return( comm_ ); }
00496 
00497 #ifdef FEI_HAVE_IOSFWD
00498    int setDbgOut(std::ostream& ostr, const char* path, const char* feiName);
00499 #else
00500    int setDbgOut(ostream& ostr, const char* path, const char* feiName);
00501 #endif
00502 
00503  private:
00504 
00505    NodeDescriptor& findNodeDescriptor(GlobalID nodeID);
00506 
00507    int writeEqn2NodeMap();
00508 
00509    int getElemNodeDescriptors(int blockIndex, int elemIndex,
00510                               NodeDescriptor** nodes);
00511 
00512    int getNodeIndices_simple(NodeDescriptor** nodes, int numNodes,
00513            int fieldID,
00514            int* scatterIndices, int& offset);
00515 
00516    int getNodeIndices_simple(NodeDescriptor** nodes, int numNodes,
00517            int fieldID,
00518            int* scatterIndices, int& offset,
00519            int* blkScatterIndices,
00520            int* blkSizes, int& blkOffset);
00521 
00522    int getNodeMajorIndices(NodeDescriptor** nodes, int numNodes,
00523                            int** fieldIDs, int* fieldsPerNode,
00524                            int* scatterIndices, int& offset);
00525 
00526    int getNodeBlkIndices(NodeDescriptor** nodes, int numNodes,
00527        int* scatterIndices, int& offset);
00528 
00529    int getNodeMajorIndices(NodeDescriptor** nodes, int numNodes,
00530                            int** fieldIDs, int* fieldsPerNode,
00531                            int* scatterIndices, int& offset,
00532          int* blkScatterIndices,
00533          int* blkSizes, int& blkOffset);
00534 
00535    int getNodeMajorIndices(NodeDescriptor** nodes, int numNodes,
00536                            std::vector<int>* fieldIDs,
00537                            std::vector<int>& fieldsPerNode,
00538                            std::vector<int>& scatterIndices);
00539 
00540    int getFieldMajorIndices(NodeDescriptor** nodes, int numNodes,
00541                             int** fieldIDs, int* fieldsPerNode,
00542                             int* scatterIndices, int& offset);
00543 
00544    int getFieldMajorIndices(NodeDescriptor** nodes, int numNodes,
00545                                           std::vector<int>* fieldIDs,
00546                                           std::vector<int>& fieldsPerNode,
00547                                           std::vector<int>& scatterIndices);
00548 
00549    void calcGlobalEqnInfo(int numLocallyOwnedNodes, 
00550         int numLocalEqns,
00551         int numLocalEqnBlks);
00552 
00553    int finalizeActiveNodes();
00554    int finalizeNodeCommMgr();
00555    bool activeNodesInitialized();
00556 
00557    int formMatrixStructure();
00558 
00559    int initElemBlockStructure();
00560    int initMultCRStructure();
00561    int initPenCRStructure();
00562    int createMatrixPosition(int row, int col, const char* callingFunction);
00563    int createMatrixPositions(int row, int numCols, int* cols,
00564            const char* callingFunction);
00565    int createMatrixPositions(fei::CSRMat& mat);
00566 
00567    int createSymmEqnStructure(std::vector<int>& scatterIndices);
00568    int createBlkSymmEqnStructure(std::vector<int>& scatterIndices);
00569 
00570    int storeElementScatterIndices(std::vector<int>& scatterIndices);
00571    int storeElementScatterIndices_noSlaves(std::vector<int>& scatterIndices);
00572    int storeElementScatterBlkIndices_noSlaves(std::vector<int>& scatterIndices);
00573 
00574    void storeLocalNodeIndices(NodeDescriptor& iNode, int iField,
00575                                      NodeDescriptor& jNode, int jField);
00576    void storeNodalColumnIndices(int eqn, NodeDescriptor& node,
00577                  int fieldID);
00578    void storeNodalRowIndices(NodeDescriptor& node, int fieldID, int eqn);
00579    void storeNodalSendIndex(NodeDescriptor& node, int fieldID, int col);
00580    void storeNodalSendIndices(NodeDescriptor& iNode, int iField,
00581                                   NodeDescriptor& jNode, int jField);
00582 
00583    int assembleReducedStructure();
00584 
00585    bool nodalEqnsAllSlaves(const NodeDescriptor* node, std::vector<int>& slaveEqns);
00586 
00587    int initializeBlkEqnMapper();
00588 
00589    int setNumNodesAndEqnsPerBlock();
00590 
00591    void destroyBlockRoster();
00592 
00593 #ifdef FEI_HAVE_IOSFWD
00594    std::ostream& dbgOut() { return( *dbgOStreamPtr_ ); }
00595 #else
00596    ostream& dbgOut() { return( *dbgOStreamPtr_ ); }
00597 #endif
00598 
00599    void addCR(int CRID,
00600        snl_fei::Constraint<GlobalID>*& cr,
00601        std::map<GlobalID,snl_fei::Constraint<GlobalID>* >& crDB);
00602 
00603    int setNodalEqnInfo();
00604    void setElemDOFEqnInfo();
00605    int setMultCREqnInfo();
00606 
00607    SNL_FEI_Structure(const SNL_FEI_Structure& src);
00608    SNL_FEI_Structure& operator=(const SNL_FEI_Structure& src);
00609 
00610    MPI_Comm comm_;
00611 
00612    int localProc_, masterProc_, numProcs_;
00613 
00614    std::map<int,int>* fieldDatabase_;
00615    std::vector<int> workarray_;
00616 
00617    std::vector<GlobalID> blockIDs_;
00618    std::vector<BlockDescriptor*> blocks_;
00619    std::vector<ConnectivityTable*> connTables_;
00620 
00621    NodeDatabase* nodeDatabase_;
00622 
00623    bool activeNodesInitialized_;
00624 
00625    std::vector<int> globalNodeOffsets_;
00626    std::vector<int> globalEqnOffsets_;
00627    std::vector<int> globalBlkEqnOffsets_;
00628 
00629    std::vector<SlaveVariable*>* slaveVars_;
00630    EqnBuffer* slaveEqns_;
00631    std::vector<int>* slvEqnNumbers_;
00632    int numSlvs_, lowestSlv_, highestSlv_;
00633    fei::FillableMat* slaveMatrix_;
00634    std::vector<int> globalNumNodesVanished_;
00635    std::vector<int> localVanishedNodeNumbers_;
00636 
00637    NodeCommMgr* nodeCommMgr_;
00638    EqnCommMgr* eqnCommMgr_;
00639    EqnCommMgr* slvCommMgr_;
00640 
00641    int numGlobalEqns_;
00642    int numLocalEqns_;
00643    int localStartRow_;
00644    int localEndRow_;
00645 
00646    int numLocalNodalEqns_;
00647    int numLocalElemDOF_;
00648    int numLocalMultCRs_;
00649 
00650    int reducedStartRow_, reducedEndRow_, numLocalReducedRows_;
00651    fei::FillableMat *Kid_, *Kdi_, *Kdd_;
00652    fei::CSRMat csrD, csrKid, csrKdi, csrKdd, tmpMat1_, tmpMat2_;
00653    int reducedEqnCounter_, reducedRHSCounter_;
00654    std::vector<int> rSlave_, cSlave_;
00655    std::vector<NodeDescriptor*> work_nodePtrs_;
00656 
00657    bool structureFinalized_;
00658    bool generateGraph_;
00659 
00660    fei::ctg_set<int>* sysMatIndices_;
00661 
00662    bool blockMatrix_;
00663    int numGlobalEqnBlks_;
00664    int numLocalEqnBlks_;
00665    int numLocalReducedEqnBlks_;
00666    int localBlkOffset_;
00667    int localReducedBlkOffset_;
00668    int globalMaxBlkSize_;
00669 
00670    int firstLocalNodeNumber_;
00671    int numGlobalNodes_;
00672 
00673    fei::ctg_set<int>* sysBlkMatIndices_;
00674    bool matIndicesDestroyed_;
00675 
00676    std::vector<int> workSpace_;
00677 
00678    snl_fei::PointBlockMap* blkEqnMapper_;
00679 
00680    std::map<GlobalID, snl_fei::Constraint<GlobalID>* > multCRs_;
00681 
00682    std::map<GlobalID, snl_fei::Constraint<GlobalID>* > penCRs_;
00683 
00684    bool checkSharedNodes_;
00685 
00686    std::string name_;
00687 
00688    int outputLevel_;
00689    bool debugOutput_;
00690    std::string dbgPath_;
00691 #ifdef FEI_HAVE_IOSFWD
00692    std::ostream* dbgOStreamPtr_;
00693 #else
00694    ostream* dbgOStreamPtr_;
00695 #endif
00696    bool setDbgOutCalled_;
00697 };
00698 
00699 #endif
00700 

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