FEI Version of the Day
fei_Lookup_Impl.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 
00044 #include <fei_macros.hpp>
00045 
00046 #include <fei_Lookup_Impl.hpp>
00047 
00048 #include <fei_VectorSpace.hpp>
00049 
00050 #include <snl_fei_Utils.hpp>
00051 #include <fei_TemplateUtils.hpp>
00052 #include <fei_CommUtils.hpp>
00053 
00054 #include <snl_fei_Constraint.hpp>
00055 
00056 #include <snl_fei_SubdMsgHandler.hpp>
00057 
00058 #undef fei_file
00059 #define fei_file "fei_Lookup_Impl.cpp"
00060 #include <fei_ErrMacros.hpp>
00061 
00062 //----------------------------------------------------------------------------
00063 fei::Lookup_Impl::Lookup_Impl(fei::SharedPtr<fei::MatrixGraph> matGraph,
00064           int nodeIDType)
00065   : matGraph_(matGraph),
00066     ptBlkMap_(NULL),
00067     vspace_(),
00068     nodeIDType_(nodeIDType),
00069     nodenumPairs_(),
00070     eqnnumPairs_(),
00071     nodenumSubdomainDB_(),
00072     databasesBuilt_(false),
00073     fieldIDs_(),
00074     fieldSizes_(),
00075     elemBlockIDs_(0, 4),
00076     fieldIDs_2D_(),
00077     workspace_()
00078 {
00079   vspace_ = matGraph_->getRowSpace();
00080   ptBlkMap_ = vspace_->getPointBlockMap();
00081 
00082   int err = buildDatabases();
00083   if (err != 0) {
00084     voidERReturn;
00085   }
00086 }
00087 
00088 //----------------------------------------------------------------------------
00089 fei::Lookup_Impl::~Lookup_Impl()
00090 {
00091   fei::destroyValues(nodenumSubdomainDB_);
00092   nodenumSubdomainDB_.clear();
00093 }
00094 
00095 //----------------------------------------------------------------------------
00096 int fei::Lookup_Impl::getEqnNumber(int nodeNumber, int fieldID)
00097 {
00098   std::map<int,fei::Record<int>*>::iterator
00099     nnp_iter = nodenumPairs_.find(nodeNumber);
00100 
00101   if (nnp_iter == nodenumPairs_.end()) return(-1);
00102 
00103   fei::Record<int>* node = (*nnp_iter).second;
00104 
00105   std::vector<int>& eqnNums = vspace_->getEqnNumbers();
00106   int* eqnNumbers = eqnNums.size() > 0 ? &eqnNums[0] : NULL;
00107   if (eqnNumbers == NULL) {
00108     throw std::runtime_error("Fatal error in fei::Lookup_Impl::getEqnNumber");
00109   }
00110   eqnNumbers += node->getOffsetIntoEqnNumbers();
00111   int offset = -1;
00112   node->getFieldMask()->getFieldEqnOffset(fieldID, offset);
00113   return(eqnNumbers[offset]);
00114 }
00115 
00116 //----------------------------------------------------------------------------
00117 int fei::Lookup_Impl::getAssociatedNodeNumber(int eqnNumber)
00118 {
00119   std::map<int,fei::Record<int>*>::iterator
00120     enp_iter = eqnnumPairs_.find(eqnNumber);
00121 
00122   if (enp_iter == eqnnumPairs_.end()) return(-1);
00123 
00124   fei::Record<int>* node = (*enp_iter).second;
00125 
00126   return( node->getNumber() );
00127 }
00128 
00129 //----------------------------------------------------------------------------
00130 int fei::Lookup_Impl::getAssociatedNodeID(int eqnNumber)
00131 {
00132   std::map<int,fei::Record<int>*>::iterator
00133     enp_iter = eqnnumPairs_.find(eqnNumber);
00134 
00135   if (enp_iter == eqnnumPairs_.end()) return(-1);
00136 
00137   fei::Record<int>* node = (*enp_iter).second;
00138 
00139   return( node->getID() );
00140 }
00141 
00142 //----------------------------------------------------------------------------
00143 int fei::Lookup_Impl::getAssociatedFieldID(int eqnNumber)
00144 {
00145   std::map<int,fei::Record<int>*>::iterator
00146     enp_iter = eqnnumPairs_.find(eqnNumber);
00147 
00148   if (enp_iter == eqnnumPairs_.end()) return(-1);
00149 
00150   fei::Record<int>* node = (*enp_iter).second;
00151 
00152   fei::FieldMask* fm = node->getFieldMask();
00153   const std::vector<int>& fieldIDs = fm->getFieldIDs();
00154   const std::vector<int>& fieldSizes = fm->getFieldSizes();
00155 
00156   const std::vector<int>& eqnNumbers = vspace_->getEqnNumbers();
00157 
00158   int baseEqnOffset = node->getOffsetIntoEqnNumbers();
00159   int numNodalEqns = fm->getNumIndices();
00160 
00161   if (baseEqnOffset + numNodalEqns > (int)eqnNumbers.size()) {
00162     throw std::runtime_error("fei::Lookup_Impl::getAssociatedFieldID ERROR, nodal eqn offset out of range.");
00163   }
00164 
00165   int offset = 0;
00166   int eqn = eqnNumbers[baseEqnOffset];
00167   while(eqn < eqnNumber && offset < numNodalEqns) {
00168     eqn = eqnNumbers[baseEqnOffset + ++offset];
00169   }
00170 
00171   if (eqn != eqnNumber) {
00172     throw std::runtime_error("fei::Lookup_Impl::getAssociatedFieldID ERROR, eqnNumber not found");
00173   }
00174 
00175   int fieldSize_total = 0;
00176   for(size_t i=0; i<fieldSizes.size(); ++i) {
00177     fieldSize_total += fieldSizes[i];
00178     if (fieldSize_total > offset) {
00179       return fieldIDs[i];
00180     }
00181   }
00182 
00183   return -1;
00184 }
00185 
00186 //----------------------------------------------------------------------------
00187 bool fei::Lookup_Impl::isInLocalElement(int nodeNumber)
00188 {
00189   std::map<int,std::vector<int>* >::iterator
00190     nns_iter = nodenumSubdomainDB_.find(nodeNumber);
00191   if (nns_iter != nodenumSubdomainDB_.end()) {
00192     return( true );
00193   }
00194 
00195   std::map<int,fei::Record<int>*>::iterator
00196     nnp_iter = nodenumPairs_.find(nodeNumber);
00197 
00198   return(nnp_iter != nodenumPairs_.end() ? true : false);
00199 }
00200 
00201 //----------------------------------------------------------------------------
00202 int fei::Lookup_Impl::getOffsetIntoBlkEqn(int blkEqn, int ptEqn)
00203 {
00204   //assume blkEqn is a node-number, for now.
00205   std::map<int,fei::Record<int>*>::iterator
00206     nnp_iter = nodenumPairs_.find(blkEqn);
00207 
00208   if (nnp_iter == nodenumPairs_.end()) return(-1);
00209 
00210   fei::Record<int>* node = (*nnp_iter).second;
00211 
00212   int eqn = vspace_->getEqnNumbers()[node->getOffsetIntoEqnNumbers()];
00213   return(ptEqn - eqn);
00214 }
00215 
00216 //----------------------------------------------------------------------------
00217 int fei::Lookup_Impl::buildDatabases()
00218 {
00219   if (databasesBuilt_) return(0);
00220 
00221   snl_fei::RecordCollection* collection = NULL;
00222   int err = vspace_->getRecordCollection(nodeIDType_, collection);
00223   if (err != 0) {
00224     //probably means that vspace_ doesn't have 'nodeIDType_', so we'll skip the
00225     //rest of this function. It's not a problem if vspace_ doesn't have nodeIDType_.
00226     return(0);
00227   }
00228 
00229   std::vector<int>& vspcEqnNumbers = vspace_->getEqnNumbers();
00230 
00231   std::vector<fei::Record<int> >& rvec = collection->getRecords();
00232 
00233   for(size_t i=0; i<rvec.size(); ++i) {
00234     fei::Record<int>* node = &rvec[i];
00235 
00236     std::pair<int,fei::Record<int>* > int_node_pair(node->getNumber(), node);
00237 
00238     nodenumPairs_.insert(int_node_pair);
00239 
00240     int numEqns = node->getFieldMask()->getNumIndices();
00241     int* eqnNumbers = &vspcEqnNumbers[0]
00242                     + node->getOffsetIntoEqnNumbers();
00243 
00244     for(int eq=0; eq<numEqns; ++eq) {
00245       std::pair<int,fei::Record<int>* > eqn_node_pair(eqnNumbers[eq], node);
00246       eqnnumPairs_.insert(eqn_node_pair);
00247     }
00248   }
00249 
00250   MPI_Comm comm = matGraph_->getRowSpace()->getCommunicator();
00251 
00252   int numLocalLagrangeConstraints = matGraph_->getLagrangeConstraints().size();
00253 
00254   int numGlobalLagrangeConstraints = 0;
00255   fei::GlobalSum(comm, numLocalLagrangeConstraints, numGlobalLagrangeConstraints);
00256 
00257   bool noconstraints = numGlobalLagrangeConstraints<1 ? true : false;
00258 
00259   fei::SharedIDs<int> subdomainIDs;
00260   fei::SharedIDs<int>& sharedIDs = vspace_->getSharedIDs_private(nodeIDType_);
00261 
00262   if (noconstraints == false) {
00263     snl_fei::SubdMsgHandler subdmsghndlr(collection, &sharedIDs, &subdomainIDs);
00264 
00265     if (vspace_->ownerPatterns_.size() > 0 && vspace_->sharerPatterns_.size() > 0) {
00266       subdmsghndlr.setSendPattern(vspace_->ownerPatterns_.find(nodeIDType_)->second);
00267       subdmsghndlr.setRecvPattern(vspace_->sharerPatterns_.find(nodeIDType_)->second);
00268       CHK_ERR( fei::exchange(comm, &subdmsghndlr) );
00269     }
00270 
00271     //Now the subdomainIDs object contains a mapping from each shared ID to a
00272     //list of processors that have that ID in their local subdomain.
00273     //So what we'll do next is run through the list of IDs in subdomainIDs and
00274     //for each ID, store the corresponding node-number in a database together
00275     //with a pointer to a list (vector) of the subdomain-processors.
00276   }
00277   else {
00278     subdomainIDs = sharedIDs;
00279   }
00280 
00281   int local_proc = fei::localProc(comm);
00282 
00283   fei::SharedIDs<int>::map_type& sdIDTable = subdomainIDs.getSharedIDs();
00284   fei::SharedIDs<int>::map_type::iterator
00285     sd_iter = sdIDTable.begin(),
00286     sd_end  = sdIDTable.end();
00287 
00288   for(int i=0; sd_iter != sd_end; ++i, ++sd_iter) {
00289     int id = sd_iter->first;
00290     std::set<int>& procList = sd_iter->second;
00291 
00292     fei::Record<int>* node = collection->getRecordWithID(id);
00293     if (node == NULL) {
00294       ERReturn(-1);
00295     }
00296 
00297     std::vector<int>* newarray = new std::vector<int>;
00298     std::set<int>::const_iterator
00299       p_iter = procList.begin(), p_end = procList.end();
00300 
00301     for(; p_iter != p_end; ++p_iter) {
00302       int proc = *p_iter;
00303       newarray->push_back(proc);
00304     }
00305 
00306     if (node->isInLocalSubdomain_) {
00307       fei::sortedListInsert(local_proc, *newarray);
00308     }
00309 
00310     nodenumSubdomainDB_.insert(std::pair<int,std::vector<int>*>(node->getNumber(), newarray));
00311   }
00312 
00313   databasesBuilt_ = true;
00314   return(0);
00315 }
00316 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends