fei_VectorSpace.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 "fei_utils.hpp"
00013 
00014 #include "fei_TemplateUtils.hpp"
00015 #include "fei_chk_mpi.hpp"
00016 #include <fei_CommUtils.hpp>
00017 #include <fei_set_shared_ids.hpp>
00018 #include "snl_fei_Utils.hpp"
00019 #include "fei_Record.hpp"
00020 #include "snl_fei_RecordCollection.hpp"
00021 #include "fei_ParameterSet.hpp"
00022 #include "snl_fei_RecordMsgHandler.hpp"
00023 #include "fei_SharedIDs.hpp"
00024 #include "fei_Pattern.hpp"
00025 #include "fei_VectorSpace.hpp"
00026 #include "fei_FieldMask.hpp"
00027 #include "snl_fei_PointBlockMap.hpp"
00028 #include "fei_LogManager.hpp"
00029 
00030 #undef fei_file
00031 #define fei_file "fei_VectorSpace.cpp"
00032 #include "fei_ErrMacros.hpp"
00033 
00034 namespace fei {
00035   class RecordAttributeCounter : public Record_Operator<int> {
00036   public:
00037     RecordAttributeCounter(int proc)
00038       : numLocalDOF_(0),
00039         numLocalIDs_(0),
00040         numLocallyOwnedIDs_(0),
00041         numRemoteSharedDOF_(0),
00042         proc_(proc)
00043     {}
00044 
00045     virtual ~RecordAttributeCounter(){}
00046 
00047     void operator()(fei::Record<int>& record)
00048     {
00049       fei::FieldMask* mask = record.getFieldMask();
00050       int owner = record.getOwnerProc();
00051 
00052       if (owner != proc_) {
00053         numRemoteSharedDOF_ += mask->getNumIndices();
00054         return;
00055       }
00056       else {
00057         ++numLocallyOwnedIDs_;
00058       }
00059 
00060       ++numLocalIDs_;
00061 
00062       int numDOF = mask->getNumIndices();
00063 
00064       numLocalDOF_ += numDOF;  
00065     }
00066 
00067     int numLocalDOF_;
00068     int numLocalIDs_;
00069     int numLocallyOwnedIDs_;
00070     int numRemoteSharedDOF_;
00071 
00072   private:
00073     int proc_;
00074   };//class Record_Operator
00075 
00076   class BlkIndexAccessor : public Record_Operator<int> {
00077   public:
00078     BlkIndexAccessor(int localProc,
00079                      int lenBlkIndices,
00080                      int* globalBlkIndices,
00081                      int* blkSizes)
00082       : numBlkIndices_(0),
00083         proc_(localProc),
00084         lenBlkIndices_(lenBlkIndices),
00085         globalBlkIndices_(globalBlkIndices),
00086         blkSizes_(blkSizes)
00087     {
00088     }
00089 
00090     BlkIndexAccessor(int lenBlkIndices,
00091                      int* globalBlkIndices,
00092                      int* blkSizes)
00093       : numBlkIndices_(0),
00094         proc_(-1),
00095         lenBlkIndices_(lenBlkIndices),
00096         globalBlkIndices_(globalBlkIndices),
00097         blkSizes_(blkSizes)
00098     {
00099     }
00100 
00101     void operator()(fei::Record<int>& record)
00102     {
00103       int owner = record.getOwnerProc();
00104       if (owner != proc_ && proc_ > -1) {
00105         return;
00106       }
00107 
00108       fei::FieldMask* mask = record.getFieldMask();
00109       int blkSize = mask->getNumIndices();
00110 
00111       if (numBlkIndices_ < lenBlkIndices_) {
00112         globalBlkIndices_[numBlkIndices_] = record.getNumber();
00113         blkSizes_[numBlkIndices_] = blkSize;
00114       }
00115 
00116       ++numBlkIndices_;
00117     }
00118 
00119     int numBlkIndices_;
00120 
00121   private:
00122     int proc_;
00123     int lenBlkIndices_;
00124     int* globalBlkIndices_;
00125     int* blkSizes_;
00126   };//class BlkIndexAccessor
00127 }//namespace snl_fei
00128 
00129 //----------------------------------------------------------------------------
00130 fei::SharedPtr<fei::VectorSpace>
00131 fei::VectorSpace::Factory::createVectorSpace(MPI_Comm comm,
00132                                              const char* name)
00133 {
00134   fei::SharedPtr<fei::VectorSpace> sptr(new fei::VectorSpace(comm, name));
00135   return(sptr);
00136 }
00137 
00138 //----------------------------------------------------------------------------
00139 fei::VectorSpace::VectorSpace(MPI_Comm comm, const char* name)
00140   : fieldMasks_(),
00141     comm_(comm),
00142     idTypes_(),
00143     fieldDatabase_(),
00144     fieldDofMap_(),
00145     recordCollections_(),
00146     sharedIDTables_(),
00147     ownerPatterns_(),
00148     sharerPatterns_(),
00149     sharedRecordsSynchronized_(true),
00150     ptBlkMap_(NULL),
00151     globalOffsets_(),
00152     globalIDOffsets_(),
00153     simpleProblem_(false),
00154     firstLocalOffset_(-1),
00155     lastLocalOffset_(-1),
00156     eqnNumbers_(),
00157     newInitData_(false),
00158     name_(),
00159     dbgprefix_("VecSpc: "),
00160     checkSharedIDs_(false)
00161 {
00162   check_version();
00163 
00164 //The initializations below are redundant with the ones above in the
00165 //initializer list. For some reason, when compiled for janus with optimization,
00166 //these pointers are not being set to NULL by the above initializers.
00167 //This was causing seg faults later when other VectorSpace methods delete these
00168 //pointers if they aren't NULL.  ABW 5/19/2004
00169 
00170   ptBlkMap_ = NULL;
00171 
00172   int numProcs = fei::numProcs(comm_);
00173   globalOffsets_.assign(numProcs+1, 0);
00174   globalIDOffsets_.assign(numProcs+1, 0);
00175 
00176   setName(name);
00177 }
00178 
00179 //----------------------------------------------------------------------------
00180 fei::VectorSpace::~VectorSpace()
00181 {
00182   int i, len = fieldMasks_.size();
00183   for(i=0; i<len; ++i) delete fieldMasks_[i];
00184 
00185   len = recordCollections_.size();
00186   for(i=0; i<len; ++i) delete recordCollections_[i];
00187 
00188   std::map<int,fei::comm_map*>::iterator
00189     o_iter = ownerPatterns_.begin(), o_end = ownerPatterns_.end();
00190   for(; o_iter != o_end; ++o_iter) {
00191     delete o_iter->second;
00192   }
00193 
00194   std::map<int,fei::comm_map*>::iterator
00195     s_iter = sharerPatterns_.begin(), s_end = sharerPatterns_.end();
00196   for(; s_iter != s_end; ++s_iter) {
00197     delete s_iter->second;
00198   }
00199 
00200   delete ptBlkMap_;
00201 }
00202 
00203 //----------------------------------------------------------------------------
00204 MPI_Comm
00205 fei::VectorSpace::getCommunicator() const
00206 {
00207   return comm_;
00208 }
00209 
00210 //----------------------------------------------------------------------------
00211 void fei::VectorSpace::setParameters(const fei::ParameterSet& paramset)
00212 {
00213   const fei::Param* param = paramset.get("name");
00214   fei::Param::ParamType ptype = param != NULL ?
00215     param->getType() : fei::Param::BAD_TYPE;
00216   if (ptype == fei::Param::STRING) {
00217     setName(param->getStringValue().c_str());
00218   }
00219 
00220   param = paramset.get("FEI_OUTPUT_LEVEL");
00221   ptype = param != NULL ? param->getType() : fei::Param::BAD_TYPE;
00222   if (ptype == fei::Param::STRING) {
00223     fei::LogManager& log_manager = fei::LogManager::getLogManager();
00224     log_manager.setOutputLevel(param->getStringValue().c_str());
00225     setOutputLevel(fei::utils::string_to_output_level(param->getStringValue()));
00226   }
00227 
00228   param = paramset.get("FEI_LOG_EQN");
00229   ptype =  param != NULL ? param->getType() : fei::Param::BAD_TYPE;  
00230   if (ptype == fei::Param::INT) {
00231     addLogEqn(param->getIntValue());
00232   }
00233 
00234   param = paramset.get("FEI_LOG_ID");
00235   ptype =  param != NULL ? param->getType() : fei::Param::BAD_TYPE;  
00236   if (ptype == fei::Param::INT) {
00237     addLogID(param->getIntValue());
00238   }
00239 
00240   param = paramset.get("FEI_CHECK_SHARED_IDS");
00241   ptype = param != NULL ? param->getType() : fei::Param::BAD_TYPE;
00242   if (ptype != fei::Param::BAD_TYPE) {
00243     if (ptype == fei::Param::BOOL) {
00244       checkSharedIDs_ = param->getBoolValue();
00245     }
00246     else if (ptype == fei::Param::INT) {
00247       checkSharedIDs_ = param->getIntValue() > 0 ? true : false;
00248     }
00249     else {
00250       checkSharedIDs_ = true;
00251     }
00252   }
00253   else {
00254     checkSharedIDs_ = false;
00255   }
00256 }
00257 
00258 //----------------------------------------------------------------------------
00259 void fei::VectorSpace::defineFields(int numFields,
00260                                     const int* fieldIDs,
00261                                     const int* fieldSizes,
00262                                     const int* fieldTypes)
00263 {
00264   if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
00265     FEI_OSTREAM& os = *output_stream_;
00266     os << dbgprefix_<<"defineFields ";
00267     for(int j=0; j<numFields; ++j) {
00268       os << "{"<<fieldIDs[j] << "," << fieldSizes[j] << "} ";
00269     }
00270     os << FEI_ENDL;
00271   }
00272 
00273   for (int i=0; i<numFields; ++i) {
00274     fieldDatabase_.insert(std::pair<int,unsigned>(fieldIDs[i], fieldSizes[i]));
00275     if (fieldIDs[i] >= 0) {
00276       if (fieldTypes != NULL) {
00277         fieldDofMap_.add_field(fieldIDs[i], fieldSizes[i], fieldTypes[i]);
00278       }
00279       else {
00280         fieldDofMap_.add_field(fieldIDs[i], fieldSizes[i]);
00281       }
00282     }
00283   }
00284 }
00285 
00286 //----------------------------------------------------------------------------
00287 void fei::VectorSpace::defineIDTypes(int numIDTypes,
00288                                      const int* idTypes)
00289 {
00290   if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
00291     FEI_OSTREAM& os = *output_stream_;
00292     os << dbgprefix_<<"defineIDTypes {";
00293     for(int j=0; j<numIDTypes; ++j) {
00294       os << idTypes[j] << " ";
00295     }
00296     os << "}"<<FEI_ENDL;
00297   }
00298 
00299   int localProc = fei::localProc(comm_);
00300   for (int i=0; i<numIDTypes; ++i) {
00301     int offset = fei::sortedListInsert(idTypes[i], idTypes_);
00302 
00303     if (offset >= 0) {
00304       recordCollections_.insert(recordCollections_.begin()+offset, new snl_fei::RecordCollection(localProc));
00305     }
00306   }
00307 }
00308 
00309 //----------------------------------------------------------------------------
00310 int fei::VectorSpace::addDOFs(int fieldID,
00311                                           int numInstancesOfThisFieldPerID,
00312                                           int idType,
00313                                           int numIDs,
00314                                           const int* IDs)
00315 {
00316   if (output_level_ > fei::BRIEF_LOGS && output_stream_ != NULL) {
00317     FEI_OSTREAM& os = *output_stream_;
00318     os << dbgprefix_<<"addDOFs, fID=" << fieldID
00319        <<", idT="<<idType <<", ninst="
00320        << numInstancesOfThisFieldPerID << " {";
00321     for(int j=0; j<numIDs; ++j) {
00322       os << IDs[j] << " ";
00323       if (j>0 && j%20==0) os << FEI_ENDL << dbgprefix_;
00324     }
00325     os << "}"<<FEI_ENDL;
00326   }
00327 
00328   if (numIDs <= 0) return(0);
00329 
00330   int idx = fei::binarySearch(idType, idTypes_);
00331   if (idx < 0) ERReturn(-1);
00332 
00333   unsigned fieldSize = getFieldSize(fieldID);
00334   recordCollections_[idx]->initRecords(fieldID, fieldSize,
00335                                        numInstancesOfThisFieldPerID,
00336                                        numIDs, IDs,
00337                                        fieldMasks_);
00338   newInitData_ = true;
00339   sharedRecordsSynchronized_ = false;
00340 
00341   return(0);
00342 }
00343 
00344 //----------------------------------------------------------------------------
00345 int fei::VectorSpace::addDOFs(int fieldID,
00346                                           int numInstancesOfThisFieldPerID,
00347                                           int idType,
00348                                           int numIDs,
00349                                           const int* IDs,
00350                                           fei::Record<int>** records)
00351 {
00352   if (output_level_ > fei::BRIEF_LOGS && output_stream_ != NULL) {
00353     FEI_OSTREAM& os = *output_stream_;
00354     os << dbgprefix_<<"addDOFs*, fID=" << fieldID
00355        <<", idT="<<idType <<", ninst="
00356        << numInstancesOfThisFieldPerID << " {";
00357     for(int j=0; j<numIDs; ++j) {
00358       os << IDs[j] << " ";
00359       if (j>0 && j%20==0) os << FEI_ENDL << dbgprefix_;
00360     }
00361     os << "}"<<FEI_ENDL;
00362   }
00363 
00364   if (numIDs <= 0) return(0);
00365 
00366   int idx = fei::binarySearch(idType, idTypes_);
00367   if (idx < 0) {
00368     FEI_OSTRINGSTREAM osstr;
00369     osstr << "fei::VectorSpace::addDOFs: error, idType " << idType
00370           << " not recognized. (idTypes need to be initialized via the"
00371           << " method VectorSpace::defineIDTypes)";
00372     throw std::runtime_error(osstr.str());
00373   }
00374 
00375   unsigned fieldSize = getFieldSize(fieldID);
00376   recordCollections_[idx]->initRecords(fieldID, fieldSize,
00377                                        numInstancesOfThisFieldPerID,
00378                                        numIDs, IDs,
00379                                        fieldMasks_,
00380                                        records);
00381   newInitData_ = true;
00382   sharedRecordsSynchronized_ = false;
00383 
00384   return(0);
00385 }
00386 
00387 //----------------------------------------------------------------------------
00388 int fei::VectorSpace::addDOFs(int idType,
00389                                           int numIDs,
00390                                           const int* IDs)
00391 {
00392   if (output_level_ > fei::BRIEF_LOGS && output_stream_ != NULL) {
00393     FEI_OSTREAM& os = *output_stream_;
00394     os << dbgprefix_<<"addDOFs idT=" <<idType<<" {";
00395     for(int j=0; j<numIDs; ++j) {
00396       os << IDs[j] << " ";
00397       if (j>0 && j%20==0) os << FEI_ENDL << dbgprefix_;
00398     }
00399     os << "}"<<FEI_ENDL;
00400   }
00401 
00402   if (numIDs <= 0) return(0);
00403 
00404   int idx = fei::binarySearch(idType, idTypes_);
00405   if (idx < 0) ERReturn(-1);
00406 
00407   recordCollections_[idx]->initRecords(numIDs, IDs, fieldMasks_);
00408 
00409   newInitData_ = true;
00410   sharedRecordsSynchronized_ = false;
00411 
00412   return(0);
00413 }
00414 
00415 //----------------------------------------------------------------------------
00416 int fei::VectorSpace::addDOFs(int idType,
00417                                           int numIDs,
00418                                           const int* IDs,
00419                                           fei::Record<int>** records)
00420 {
00421   if (output_level_ > fei::BRIEF_LOGS && output_stream_ != NULL) {
00422     FEI_OSTREAM& os = *output_stream_;
00423     os << dbgprefix_<<"addDOFs* idT=" <<idType<<" {";
00424     for(int j=0; j<numIDs; ++j) {
00425       os << IDs[j] << " ";
00426       if (j>0 && j%20==0) os << FEI_ENDL << dbgprefix_;
00427     }
00428     os << "}"<<FEI_ENDL;
00429   }
00430 
00431   if (numIDs <= 0) return(0);
00432 
00433   int idx = fei::binarySearch(idType, idTypes_);
00434   if (idx < 0) ERReturn(-1);
00435 
00436   recordCollections_[idx]->initRecords(numIDs, IDs,
00437                                        fieldMasks_, records);
00438 
00439   newInitData_ = true;
00440   sharedRecordsSynchronized_ = false;
00441 
00442   return(0);
00443 }
00444 
00445 //----------------------------------------------------------------------------
00446 int fei::VectorSpace::initSharedIDs(int numShared,
00447                                     int idType,
00448                                     const int* sharedIDs,
00449                                     const int* numSharingProcsPerID,
00450                                     const int* sharingProcs)
00451 {
00452   if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
00453     FEI_OSTREAM& os = *output_stream_;
00454     os << dbgprefix_<<"initSharedIDs n=" <<numShared<<", idT="<<idType<< FEI_ENDL;
00455     int offset = 0;
00456     for(int ns=0; ns<numShared; ++ns) {
00457       os << dbgprefix_<<"#sharedID="<<sharedIDs[ns] << ", nprocs=" << numSharingProcsPerID[ns] << ", procs: ";
00458       for(int sp=0; sp<numSharingProcsPerID[ns]; ++sp) {
00459         os << sharingProcs[offset++] << " ";
00460       }
00461       os << FEI_ENDL;
00462     }
00463     os << FEI_ENDL;
00464   }
00465 
00466   if (numShared == 0) return(0);
00467 
00468   int idx = fei::binarySearch(idType, idTypes_);
00469   if (idx < 0) ERReturn(-1);
00470 
00471   fei::SharedIDs<int>& shIDs = getSharedIDs_private(idType);
00472 
00473   int offset = 0;
00474   for(int i=0; i<numShared; ++i) {
00475     shIDs.addSharedID(sharedIDs[i], numSharingProcsPerID[i],
00476                                 &(sharingProcs[offset]));
00477     offset += numSharingProcsPerID[i];
00478 
00479     fei::Record<int>* rec = recordCollections_[idx]->getRecordWithID(sharedIDs[i]);
00480     if (rec == NULL) {
00481       CHK_ERR( addDOFs(idType, 1, &(sharedIDs[i])) );
00482     }
00483   }
00484 
00485   newInitData_ = true;
00486   sharedRecordsSynchronized_ = false;
00487 
00488   return(0);
00489 }
00490 
00491 //----------------------------------------------------------------------------
00492 int fei::VectorSpace::initSharedIDs(int numShared,
00493                                     int idType,
00494                                     const int* sharedIDs,
00495                                     const int* numSharingProcsPerID,
00496                                     const int* const* sharingProcs)
00497 {
00498   if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
00499     FEI_OSTREAM& os = *output_stream_;
00500     os << dbgprefix_<<"initSharedIDs n=" <<numShared<<", idT="<<idType<< FEI_ENDL;
00501     for(int ns=0; ns<numShared; ++ns) {
00502       os << dbgprefix_<<"#sharedID="<<sharedIDs[ns] << ", nprocs=" << numSharingProcsPerID[ns] << ", procs: ";
00503       for(int sp=0; sp<numSharingProcsPerID[ns]; ++sp) {
00504         os << sharingProcs[ns][sp] << " ";
00505       }
00506       os << FEI_ENDL;
00507     }
00508     os << FEI_ENDL;
00509   }
00510 
00511   if (numShared == 0) return(0);
00512 
00513   fei::SharedIDs<int>& shIDs = getSharedIDs_private(idType);
00514 
00515   int idx = fei::binarySearch(idType, idTypes_);
00516   if (idx < 0) ERReturn(-1);
00517 
00518   for(int i=0; i<numShared; ++i) {
00519     shIDs.addSharedID(sharedIDs[i], numSharingProcsPerID[i],
00520                                  sharingProcs[i]);
00521 
00522     fei::Record<int>* rec = recordCollections_[idx]->getRecordWithID(sharedIDs[i]);
00523     if (rec == NULL) {
00524       CHK_ERR( addDOFs(idType, 1, &(sharedIDs[i])) );
00525     }
00526   }
00527 
00528   newInitData_ = true;
00529   sharedRecordsSynchronized_ = false;
00530 
00531   return(0);
00532 }
00533 
00534 //----------------------------------------------------------------------------
00535 int fei::VectorSpace::addVectorSpace(fei::VectorSpace* inputSpace)
00536 {
00537   idTypes_ = inputSpace->idTypes_;
00538 
00539   std::map<int,unsigned>::const_iterator
00540     f_iter = inputSpace->fieldDatabase_.begin(),
00541     f_end = inputSpace->fieldDatabase_.end();
00542 
00543   for(; f_iter != f_end; ++f_iter) {
00544     const std::pair<const int,unsigned>& fpair = *f_iter;
00545     int fieldsize = fpair.second;
00546     defineFields(1, &(fpair.first), &fieldsize);
00547   }
00548 
00549   int i, len = inputSpace->fieldMasks_.size();
00550   fieldMasks_.resize(len);
00551   for(i=0; i<len; ++i) {
00552     fieldMasks_[i] = new fei::FieldMask(*(inputSpace->fieldMasks_[i]));
00553   }
00554 
00555   len = inputSpace->recordCollections_.size();
00556   recordCollections_.resize(len);
00557   for(i=0; i<len; ++i) {
00558     recordCollections_[i] =
00559       new snl_fei::RecordCollection(*(inputSpace->recordCollections_[i]));
00560   }
00561 
00562   sharedIDTables_ = inputSpace->sharedIDTables_;
00563 
00564   newInitData_ = true;
00565   sharedRecordsSynchronized_ = false;
00566 
00567   return(0);
00568 }
00569 
00570 //----------------------------------------------------------------------------
00571 fei::SharedIDs<int>& fei::VectorSpace::getSharedIDs_private(int idType)
00572 {
00573   std::map<int,fei::SharedIDs<int> >::iterator
00574     iter = sharedIDTables_.find(idType);
00575 
00576   if (iter == sharedIDTables_.end()) {
00577     iter = sharedIDTables_.insert(std::make_pair(idType,fei::SharedIDs<int>())).first;
00578   }
00579 
00580   return iter->second;
00581 }
00582 
00583 //----------------------------------------------------------------------------
00584 void fei::VectorSpace::compute_shared_ids()
00585 {
00586   //For each ID-type, call the fei::set_shared_ids function, which
00587   //takes a RecordCollection as input and fills a SharedIDs object
00588   //with mappings for which processors share each ID that appears
00589   //on multiple processors.
00590   for(size_t i=0; i<idTypes_.size(); ++i) {
00591     snl_fei::RecordCollection* records = recordCollections_[i];
00592     fei::SharedIDs<int>& shIDs = getSharedIDs_private(idTypes_[i]);
00593 
00594     fei::set_shared_ids(comm_, *records, shIDs);
00595   }
00596 }
00597 
00598 //----------------------------------------------------------------------------
00599 int fei::VectorSpace::initComplete()
00600 {
00601   if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
00602     FEI_OSTREAM& os = *output_stream_;
00603     os <<dbgprefix_<< "initComplete" << FEI_ENDL;
00604   }
00605 
00606   simpleProblem_ = (fieldMasks_.size()==1) && (fieldMasks_[0]->getNumFields()==1);
00607 
00608   //we need to know if any processor has newInitData_.
00609   int localInitData = newInitData_ ? 1 : 0;
00610   int globalInitData = 0;
00611   CHK_ERR( fei::GlobalMax(comm_, localInitData, globalInitData) );
00612   newInitData_ = globalInitData > 0 ? true : false;
00613 
00614   if (!newInitData_) {
00615     return(0);
00616   }
00617 
00618   //if running on multiple procs and initSharedIDs(...) has not been
00619   //called, then we need to specifically figure out which IDs are
00620   //shared and which procs share them.
00621 
00622   bool need_to_compute_shared_ids = false;
00623   if (fei::numProcs(comm_) > 1 && sharedIDTables_.size() < 1) {
00624     need_to_compute_shared_ids = true;
00625   }
00626 
00627   if (need_to_compute_shared_ids) {
00628     compute_shared_ids();
00629   }
00630 
00631   //setOwners_lowestSharing is a local operation (no communication), which
00632   //assumes that each processor holds CORRECT (globally symmetric)
00633   //shared-id/sharing-proc tables. No correctness-checking is performed here.
00634   setOwners_lowestSharing();
00635 
00636   //synchronizeSharedRecords ensures that each sharing processor has the same
00637   //view of the shared records, with respect to the layout of fields, which
00638   //determines how many DOFs and equation-numbers reside at each ID.
00639   //This involves inter-processor communication.
00640   if ( synchronizeSharedRecords() != 0) {
00641     return(-1);
00642   }
00643 
00644   //calculateGlobalIndices is also a global operation.
00645   if ( calculateGlobalIndices() != 0) {
00646     return(-1);
00647   }
00648 
00649   //finally we need to exchange global indices for shared records. i.e., the
00650   //processors that own shared records, need to send the global indices for
00651   //those records to the sharing-but-not-owning processors.
00652   if (fei::numProcs(comm_) > 1) {
00653     if ( exchangeGlobalIndices() != 0) {
00654       return(-1);
00655     }
00656   }
00657 
00658   newInitData_ = false;
00659 
00660   return(0);
00661 }
00662 
00663 //----------------------------------------------------------------------------
00664 int fei::VectorSpace::getGlobalIndex(int idType,
00665                                      int ID,
00666                                      int fieldID,
00667                                      int fieldOffset,
00668                                      int whichComponentOfField,
00669                                      int& globalIndex)
00670 {
00671   int idindex = fei::binarySearch(idType, idTypes_);
00672   if (idindex < 0) return(-1);
00673 
00674   unsigned fieldSize = 0;
00675   if (fieldOffset > 0) {
00676     fieldSize = getFieldSize(fieldID);
00677   }
00678 
00679   try {
00680     globalIndex = recordCollections_[idindex]->getGlobalIndex(ID,
00681                                                          fieldID,
00682                                                          fieldSize,
00683                                                          fieldOffset,
00684                                                          whichComponentOfField,
00685                                                          &eqnNumbers_[0]);
00686   }
00687   catch (std::runtime_error& exc) {
00688     FEI_OSTRINGSTREAM osstr;
00689     osstr << "VectorSpace::getGlobalIndex caught exception: " << exc.what();
00690     FEI_CERR << osstr.str()<<FEI_ENDL;
00691     ERReturn(-1);
00692   }
00693 
00694   return(0);
00695 }
00696 
00697 //----------------------------------------------------------------------------
00698 int fei::VectorSpace::getGlobalIndex(int idType,
00699                                      int ID,
00700                                      int fieldID,
00701                                      int& globalIndex)
00702 {
00703   return( getGlobalIndex(idType, ID, fieldID, 0, 0, globalIndex) );
00704 }
00705 
00706 //----------------------------------------------------------------------------
00707 int fei::VectorSpace::getGlobalBlkIndex(int idType,
00708                                         int ID,
00709                                         int& globalBlkIndex)
00710 {
00711   int idindex = fei::binarySearch(idType, idTypes_);
00712   if (idindex < 0) return(-1);
00713 
00714   CHK_ERR(recordCollections_[idindex]->getGlobalBlkIndex(ID, globalBlkIndex));
00715 
00716   return(0);
00717 }
00718 
00719 //----------------------------------------------------------------------------
00720 int fei::VectorSpace::getGlobalIndices(int numIDs,
00721                                        const int* IDs,
00722                                        int idType,
00723                                        int fieldID,
00724                                        int* globalIndices)
00725 {
00726   int idindex = fei::binarySearch(idType, idTypes_);
00727   if (idindex < 0) return(-1);
00728 
00729   unsigned fieldSize = getFieldSize(fieldID);
00730   int offset = 0;
00731 
00732   for(int i=0; i<numIDs; ++i) {
00733     try {
00734       globalIndices[offset] =
00735         recordCollections_[idindex]->getGlobalIndex(IDs[i],
00736                                                     fieldID,
00737                                                     fieldSize,
00738                                                     0, 0,
00739                                                     &eqnNumbers_[0]);
00740       if (fieldSize>1) {
00741         int eqn = globalIndices[offset];
00742         for(unsigned j=1; j<fieldSize; ++j) {
00743           globalIndices[offset+j] = eqn+j;
00744         }
00745       }
00746     }
00747     catch (...) {
00748       for(unsigned j=0; j<fieldSize; ++j) {
00749         globalIndices[offset+j] = -1;
00750       }
00751     }
00752 
00753     offset += fieldSize;
00754   }
00755 
00756   return(0);
00757 }
00758 
00759 //----------------------------------------------------------------------------
00760 int fei::VectorSpace::getGlobalBlkIndices(int numIDs,
00761                                           const int* IDs,
00762                                           int idType,
00763                                           int* globalBlkIndices)
00764 {
00765   int err;
00766   int idindex = fei::binarySearch(idType, idTypes_);
00767   if (idindex < 0) return(-1);
00768 
00769   int offset = 0;
00770 
00771   for(int i=0; i<numIDs; ++i) {
00772     err = recordCollections_[idindex]->getGlobalBlkIndex(IDs[i],
00773                                                      globalBlkIndices[offset]);
00774     if (err != 0) {
00775       globalBlkIndices[offset] = -1;
00776     }
00777     ++offset;
00778   }
00779 
00780   return(0);
00781 }
00782 
00783 //----------------------------------------------------------------------------
00784 int fei::VectorSpace::getGlobalIndices(int numIDs,
00785                                        const int* IDs,
00786                                        const int* idTypes,
00787                                        const int* fieldIDs,
00788                                        int* globalIndices)
00789 {
00790   int err;
00791   int offset = 0;
00792   for(int i=0; i<numIDs; ++i) {
00793     unsigned fieldSize = getFieldSize(fieldIDs[i]);
00794     err = getGlobalIndex(idTypes[i], IDs[i], fieldIDs[i], 0, 0,
00795                          globalIndices[offset]);
00796     if (err) {
00797       for(unsigned j=1; j<fieldSize; ++j) {
00798         globalIndices[offset+j] = -1;
00799       }
00800     }
00801     else {
00802       if (fieldSize>1) {
00803         int eqn = globalIndices[offset];
00804         for(unsigned j=1; j<fieldSize; ++j) {
00805           globalIndices[offset+j] = eqn+j;
00806         }
00807       }
00808     }
00809     offset += fieldSize;
00810   }
00811 
00812   return(0);
00813 }
00814 
00815 //----------------------------------------------------------------------------
00816 void fei::VectorSpace::getGlobalBlkIndices(const fei::Pattern* pattern,
00817                                            const fei::Record<int>*const* records,
00818                                            std::vector<int>& indices)
00819 {
00820   int numRecords = pattern->getNumIDs();
00821   indices.resize(numRecords);
00822   int numIndices;
00823   getGlobalBlkIndices(numRecords, records, numRecords, &indices[0],
00824                       numIndices);
00825 }
00826 
00827 //----------------------------------------------------------------------------
00828 void fei::VectorSpace::getGlobalIndices(const fei::Pattern* pattern,
00829                                         const fei::Record<int>*const* records,
00830                                         std::vector<int>& indices)
00831 {
00832   int numRecords = pattern->getNumIDs();
00833   int numIndices = pattern->getNumIndices();
00834   indices.resize(numIndices);
00835   int* indices_ptr = &indices[0];
00836 
00837   fei::Pattern::PatternType pType = pattern->getPatternType();
00838 
00839   if (pType == fei::Pattern::GENERAL ||
00840       pType == fei::Pattern::SINGLE_IDTYPE) {
00841     const int* numFieldsPerID = pattern->getNumFieldsPerID();
00842     const int* fieldIDs = pattern->getFieldIDs();
00843     int totalNumFields = pattern->getTotalNumFields();
00844 
00845     std::vector<int> fieldSizes(totalNumFields);
00846 
00847     for(int j=0; j<totalNumFields; ++j) {
00848       fieldSizes[j] = getFieldSize(fieldIDs[j]);
00849     }
00850 
00851     getGlobalIndices(numRecords, records, numFieldsPerID,
00852                      fieldIDs, &(fieldSizes[0]),
00853                      numIndices, indices_ptr, numIndices);
00854   }
00855   else if (pType == fei::Pattern::SIMPLE) {
00856     const int* fieldIDs = pattern->getFieldIDs();
00857 
00858     int fieldID = fieldIDs[0];
00859     unsigned fieldSize = getFieldSize(fieldID);
00860 
00861     getGlobalIndices(numRecords, records,
00862                      fieldID, fieldSize,
00863                      numIndices, indices_ptr, numIndices);
00864   }
00865   else if (pType == fei::Pattern::NO_FIELD) {
00866     getGlobalBlkIndices(numRecords, records, numIndices, indices_ptr, numIndices);
00867   }
00868 }
00869 
00870 //----------------------------------------------------------------------------
00871 void fei::VectorSpace::getGlobalIndices(int numRecords,
00872                                         const fei::Record<int>*const* records,
00873                                         int fieldID,
00874                                         int fieldSize,
00875                                         int indicesAllocLen,
00876                                         int* indices,
00877                                         int& numIndices)
00878 {
00879   numIndices = 0;
00880   int* eqnPtr = &eqnNumbers_[0];
00881 
00882   int len = numRecords;
00883   if (len*fieldSize >= indicesAllocLen) {
00884     len = indicesAllocLen/fieldSize;
00885   }
00886 
00887   if (fieldSize == 1 && simpleProblem_) {
00888     for(int i=0; i<len; ++i) {
00889       const fei::Record<int>* record = records[i];
00890       indices[numIndices++] = *(eqnPtr+record->getOffsetIntoEqnNumbers());
00891     }    
00892     return;
00893   }
00894 
00895   if (fieldSize == 1) {
00896     for(int i=0; i<len; ++i) {
00897       const fei::Record<int>* record = records[i];
00898 
00899       int eqnOffset = 0, numInstances = 0;
00900       int* eqnNumbers = eqnPtr+record->getOffsetIntoEqnNumbers();
00901 
00902       const fei::FieldMask* fieldMask = record->getFieldMask();
00903       fieldMask->getFieldEqnOffset(fieldID, eqnOffset, numInstances);
00904       indices[numIndices++] = eqnNumbers[eqnOffset];
00905     }
00906   }
00907   else {
00908     for(int i=0; i<len; ++i) {
00909       const fei::Record<int>* record = records[i];
00910 
00911       int* eqnNumbers = eqnPtr+record->getOffsetIntoEqnNumbers();
00912 
00913       int eqnOffset = 0, numInstances = 0;
00914       if (!simpleProblem_) {
00915         const fei::FieldMask* fieldMask = record->getFieldMask();
00916         fieldMask->getFieldEqnOffset(fieldID, eqnOffset, numInstances);
00917       }
00918 
00919       for(int fs=0; fs<fieldSize; ++fs) {
00920         indices[numIndices++] = eqnNumbers[eqnOffset+fs];
00921       }
00922     }
00923   }
00924 }
00925 
00926 //----------------------------------------------------------------------------
00927 void fei::VectorSpace::getGlobalIndices(int numRecords,
00928                                         const fei::Record<int>*const* records,
00929                                         const int* numFieldsPerID,
00930                                         const int* fieldIDs,
00931                                         const int* fieldSizes,
00932                                         int indicesAllocLen,
00933                                         int* indices,
00934                                         int& numIndices)
00935 {
00936   numIndices = 0;
00937   int fld_offset = 0;
00938   int* eqnPtr = &eqnNumbers_[0];
00939 
00940   for(int i=0; i<numRecords; ++i) {
00941     const fei::Record<int>* record = records[i];
00942 
00943     const fei::FieldMask* fieldMask = record->getFieldMask();
00944     int* eqnNumbers = eqnPtr + record->getOffsetIntoEqnNumbers();
00945 
00946     for(int nf=0; nf<numFieldsPerID[i]; ++nf) {
00947       int eqnOffset = 0, numInstances = 0;
00948       if (!simpleProblem_) {
00949         fieldMask->getFieldEqnOffset(fieldIDs[fld_offset], eqnOffset, numInstances);
00950       }
00951 
00952       for(int fs=0; fs<fieldSizes[fld_offset]; ++fs) {
00953         indices[numIndices++] = eqnNumbers[eqnOffset+fs];
00954       }
00955 
00956       ++fld_offset;
00957     }
00958   }
00959 }
00960 
00961 //----------------------------------------------------------------------------
00962 void fei::VectorSpace::getGlobalBlkIndices(int numRecords,
00963                                            const fei::Record<int>*const* records,
00964                                            int indicesAllocLen,
00965                                            int* indices,
00966                                            int& numIndices)
00967 {
00968   numIndices = 0;
00969   for(int i=0; i<numRecords; ++i) {
00970     if (numIndices < indicesAllocLen) {
00971       indices[numIndices++] = records[i]->getNumber();
00972     }
00973     else ++numIndices;
00974   }
00975 }
00976 
00977 //----------------------------------------------------------------------------
00978 int fei::VectorSpace::getGlobalIndex(int idType,
00979                                      int ID,
00980                                      int& globalIndex)
00981 {
00982   int idindex = fei::binarySearch(idType, idTypes_);
00983   if (idindex < 0) return(-1);
00984 
00985   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
00986   if (record == NULL) {
00987     ERReturn(-1);
00988   }
00989 
00990   const int* eqnNums = &eqnNumbers_[0]
00991                      + record->getOffsetIntoEqnNumbers();
00992 
00993   if (eqnNums != NULL) { globalIndex = eqnNums[0]; return(0); }
00994   return(-1);
00995 }
00996 
00997 //----------------------------------------------------------------------------
00998 int fei::VectorSpace::getNumDegreesOfFreedom(int idType,
00999                                              int ID)
01000 {
01001   int idindex = fei::binarySearch(idType, idTypes_);
01002   if (idindex < 0) return(0);
01003 
01004   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
01005   if (record == NULL) {
01006     ERReturn(-1);
01007   }
01008 
01009   return( record->getFieldMask()->getNumIndices() );
01010 }
01011 
01012 //----------------------------------------------------------------------------
01013 int fei::VectorSpace::getNumFields()
01014 {
01015   return(fieldDatabase_.size());
01016 }
01017 
01018 //----------------------------------------------------------------------------
01019 void fei::VectorSpace::getFields(std::vector<int>& fieldIDs)
01020 {
01021   unsigned numFields = fieldDatabase_.size();
01022 
01023   fieldIDs.resize(numFields);
01024 
01025   fei::copyKeysToArray<std::map<int,unsigned> >(fieldDatabase_, numFields,
01026                                                     &fieldIDs[0]);
01027 }
01028 
01029 //----------------------------------------------------------------------------
01030 size_t fei::VectorSpace::getNumIDTypes()
01031 {
01032   return(idTypes_.size());
01033 }
01034 
01035 //----------------------------------------------------------------------------
01036 void fei::VectorSpace::getIDTypes(std::vector<int>& idTypes) const
01037 {
01038   size_t numIDTypes = idTypes_.size();
01039   idTypes.resize(numIDTypes);
01040   for(size_t i=0; i<numIDTypes; ++i) idTypes[i] = idTypes_[i];
01041 }
01042 
01043 //----------------------------------------------------------------------------
01044 int fei::VectorSpace::getNumFields(int idType,
01045                                    int ID)
01046 {
01047   int idindex = fei::binarySearch(idType, idTypes_);
01048   if (idindex < 0) return(0);
01049 
01050   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
01051   if (record == NULL) {
01052     ERReturn(-1);
01053   }
01054 
01055   return( record->getFieldMask()->getNumFields() );
01056 }
01057 
01058 //----------------------------------------------------------------------------
01059 bool fei::VectorSpace::isLocal(int idType, int ID)
01060 {
01061   int idindex = fei::binarySearch(idType, idTypes_);
01062   if (idindex < 0) return(false);
01063 
01064   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
01065   if (record == NULL) {
01066     return(false);
01067   }
01068 
01069   return(true);
01070 }
01071 
01072 //----------------------------------------------------------------------------
01073 bool fei::VectorSpace::isLocallyOwned(int idType, int ID)
01074 {
01075   int idindex = fei::binarySearch(idType, idTypes_);
01076   if (idindex < 0) return(false);
01077 
01078   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
01079   if (record == NULL) {
01080     return(false);
01081   }
01082   if (record->getOwnerProc() == fei::localProc(comm_)) {
01083     return(true);
01084   }
01085 
01086   return(false);
01087 }
01088 
01089 //----------------------------------------------------------------------------
01090 unsigned fei::VectorSpace::getFieldSize(int fieldID)
01091 {
01092   std::map<int,unsigned>::const_iterator
01093     f_iter = fieldDatabase_.find(fieldID);
01094 
01095   if (f_iter == fieldDatabase_.end()) {
01096     FEI_OSTRINGSTREAM osstr;
01097     osstr << "fei::VectorSpace";
01098     if (name_.length() > 0) {
01099       osstr << "(name: "<<name_<<")";
01100     }
01101     osstr << "::getFieldSize: fieldID " << fieldID << " not found.";
01102     throw std::runtime_error(osstr.str());
01103   }
01104 
01105   return((*f_iter).second);
01106 }
01107 
01108 //----------------------------------------------------------------------------
01109 void fei::VectorSpace::getFields(int idType, int ID, std::vector<int>& fieldIDs)
01110 {
01111   int idindex = fei::binarySearch(idType, idTypes_);
01112   if (idindex < 0) {
01113     fieldIDs.resize(0);
01114     return;
01115   }
01116 
01117   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
01118   if (record == NULL) {
01119     voidERReturn;
01120   }
01121 
01122   fei::FieldMask* fieldMask = record->getFieldMask();
01123   std::vector<int>& maskFieldIDs = fieldMask->getFieldIDs();
01124   int numFields = maskFieldIDs.size();
01125   fieldIDs.resize(numFields);
01126   for(int i=0; i<numFields; ++i) {
01127     fieldIDs[i] = maskFieldIDs[i];
01128   }
01129 }
01130 
01131 //----------------------------------------------------------------------------
01132 void fei::VectorSpace::getGlobalIndexOffsets(std::vector<int>& globalOffsets) const
01133 {
01134   globalOffsets = globalOffsets_;
01135 }
01136 
01137 //----------------------------------------------------------------------------
01138 void fei::VectorSpace::getGlobalBlkIndexOffsets(std::vector<int>& globalBlkOffsets) const
01139 {
01140   globalBlkOffsets = globalIDOffsets_;
01141 }
01142 
01143 //----------------------------------------------------------------------------
01144 int fei::VectorSpace::getOwnerProcPtIndex(int globalIndex)
01145 {
01146   if (globalIndex < 0) return(-1);
01147 
01148   unsigned len = globalOffsets_.size();
01149   for(int i=0; i<(int)(len-1); ++i) {
01150     if (globalIndex < globalOffsets_[i+1]) {
01151       return(i);
01152     }
01153   }
01154 
01155   return(-1);
01156 }
01157 
01158 //----------------------------------------------------------------------------
01159 int fei::VectorSpace::getOwnerProcBlkIndex(int globalIndex)
01160 {
01161   if (globalIndex < 0) return(-1);
01162 
01163   unsigned len = globalOffsets_.size();
01164   for(int i=0; i<(int)(len-1); ++i) {
01165     if (globalIndex < globalIDOffsets_[i+1]) {
01166       return(i);
01167     }
01168   }
01169 
01170   return(-1);
01171 }
01172 
01173 //----------------------------------------------------------------------------
01174 int fei::VectorSpace::getNumOwnedAndSharedIDs(int idType)
01175 {
01176   int idx = fei::binarySearch(idType, idTypes_);
01177   if (idx < 0) return(0);
01178 
01179   return( recordCollections_[idx]->getNumRecords() );
01180 }
01181 
01182 //----------------------------------------------------------------------------
01183 int fei::VectorSpace::getNumOwnedIDs(int idType)
01184 {
01185   int idx = fei::binarySearch(idType, idTypes_);
01186   if (idx < 0) return(0);
01187 
01188   fei::RecordAttributeCounter attrCounter(fei::localProc(comm_));
01189   runRecords(attrCounter);
01190 
01191   return( attrCounter.numLocallyOwnedIDs_ );
01192 }
01193 
01194 //----------------------------------------------------------------------------
01195 int fei::VectorSpace::getNumSharedIDs(int idType, int& numShared)
01196 {
01197   numShared = sharedIDTables_[idType].getSharedIDs().size();
01198 
01199   return(0);
01200 }
01201 
01202 //----------------------------------------------------------------------------
01203 int fei::VectorSpace::getOwnedAndSharedIDs(int idType,
01204                                       int lenList,
01205                                       int* IDs,
01206                                       int& numLocalIDs)
01207 {
01208   int idx = fei::binarySearch(idType, idTypes_);
01209   if (idx < 0) return(-1);
01210 
01211   snl_fei::RecordCollection* records = recordCollections_[idx];
01212 
01213   snl_fei::RecordCollection::map_type& rmap = records->getRecords();
01214 
01215   numLocalIDs = rmap.size();
01216 
01217   fei::copyKeysToArray(rmap, lenList, IDs);
01218 
01219   return(0);
01220 }
01221 
01222 //----------------------------------------------------------------------------
01223 int fei::VectorSpace::getOwnedIDs(int idType,
01224                                          int lenList,
01225                                          int* IDs,
01226                                          int& numLocalIDs)
01227 {
01228   int idx = fei::binarySearch(idType, idTypes_);
01229   if (idx < 0) return(-1);
01230 
01231   snl_fei::RecordCollection* records = recordCollections_[idx];
01232 
01233   snl_fei::RecordCollection::map_type& rmap = records->getRecords();
01234 
01235   numLocalIDs = 0;
01236 
01237   snl_fei::RecordCollection::map_type::iterator
01238     r_iter = rmap.begin(),
01239     r_end  = rmap.end();
01240 
01241   for(; r_iter != r_end; ++r_iter) {
01242     fei::Record<int>& thisrecord = *((*r_iter).second);
01243 
01244     if (thisrecord.getOwnerProc() == fei::localProc(comm_)) {
01245       if (numLocalIDs < lenList) {
01246         IDs[numLocalIDs] = thisrecord.getID();
01247       }
01248       ++numLocalIDs;
01249     }
01250   }
01251 
01252   return(0);
01253 }
01254 
01255 //----------------------------------------------------------------------------
01256 int fei::VectorSpace::getNumIndices_SharedAndOwned() const
01257 {
01258   return(eqnNumbers_.size());
01259 }
01260 
01261 //----------------------------------------------------------------------------
01262 int fei::VectorSpace::getIndices_SharedAndOwned(std::vector<int>& globalIndices) const
01263 {
01264   if (eqnNumbers_.size() == 0) {
01265     globalIndices.resize(0);
01266     return(0);
01267   }
01268 
01269   size_t numIndices = eqnNumbers_.size();
01270   const int* indicesPtr = &eqnNumbers_[0];
01271 
01272   globalIndices.resize(numIndices);
01273   for(size_t i=0; i<numIndices; ++i) {
01274     globalIndices[i] = indicesPtr[i];
01275   }
01276 
01277   return(0);
01278 }
01279 
01280 //----------------------------------------------------------------------------
01281 int fei::VectorSpace::getNumBlkIndices_SharedAndOwned(int& numBlkIndices) const
01282 {
01283   numBlkIndices = 0;
01284   for(size_t i=0; i<recordCollections_.size(); ++i) {
01285     numBlkIndices += recordCollections_[i]->getNumRecords();
01286   }
01287 
01288   return(0);
01289 }
01290 
01291 //----------------------------------------------------------------------------
01292 int fei::VectorSpace::getBlkIndices_SharedAndOwned(int lenBlkIndices,
01293                                                    int* globalBlkIndices,
01294                                                    int* blkSizes,
01295                                                    int& numBlkIndices)
01296 {
01297   if (!sharedRecordsSynchronized_) {
01298     numBlkIndices = 0;
01299     return(-1);
01300   }
01301 
01302   fei::BlkIndexAccessor blkIndAccessor(lenBlkIndices,
01303                                   globalBlkIndices, blkSizes);
01304   runRecords(blkIndAccessor);
01305 
01306   numBlkIndices = blkIndAccessor.numBlkIndices_;
01307 
01308   return(0);
01309 }
01310 
01311 //----------------------------------------------------------------------------
01312 int fei::VectorSpace::getGlobalNumIndices() const
01313 {
01314   if (globalOffsets_.size() < 1) return(0);
01315   return(globalOffsets_[globalOffsets_.size()-1]);
01316 }
01317 
01318 //----------------------------------------------------------------------------
01319 int fei::VectorSpace::getNumIndices_Owned() const
01320 {
01321   if (!sharedRecordsSynchronized_) {
01322     return(-1);
01323   }
01324 
01325   int localProc = fei::localProc(comm_);
01326   int numIndices = globalOffsets_[localProc+1]-globalOffsets_[localProc];
01327 
01328   return(numIndices);
01329 }
01330 
01331 //----------------------------------------------------------------------------
01332 int fei::VectorSpace::getIndices_Owned(std::vector<int>& globalIndices) const
01333 {
01334   if (!sharedRecordsSynchronized_) {
01335     globalIndices.clear();
01336     return(-1);
01337   }
01338 
01339   int localProc = fei::localProc(comm_);
01340   int numIndices = globalOffsets_[localProc+1]-globalOffsets_[localProc];
01341 
01342   globalIndices.resize(numIndices);
01343 
01344   int firstOffset = globalOffsets_[localProc];
01345   for(int i=0; i<numIndices; ++i) {
01346     globalIndices[i] = firstOffset+i;
01347   }
01348 
01349   return(0);
01350 }
01351 
01352 //----------------------------------------------------------------------------
01353 int fei::VectorSpace::getIndices_Owned(int lenIndices, int* globalIndices, int& numIndices) const
01354 {
01355   if (!sharedRecordsSynchronized_) {
01356     numIndices = 0;
01357     return(-1);
01358   }
01359 
01360   int localProc = fei::localProc(comm_);
01361   numIndices = globalOffsets_[localProc+1]-globalOffsets_[localProc];
01362 
01363   int len = lenIndices >= numIndices ? numIndices : lenIndices;
01364 
01365   int firstOffset = globalOffsets_[localProc];
01366   for(int i=0; i<len; ++i) {
01367     globalIndices[i] = firstOffset+i;
01368   }
01369 
01370   return(0);
01371 }
01372 
01373 //----------------------------------------------------------------------------
01374 int fei::VectorSpace::getNumBlkIndices_Owned() const
01375 {
01376   if (!sharedRecordsSynchronized_) {
01377     return(-1);
01378   }
01379 
01380   int localProc = fei::localProc(comm_);
01381   int numBlkIndices = globalIDOffsets_[localProc+1]-globalIDOffsets_[localProc];
01382 
01383   return(numBlkIndices);
01384 }
01385 
01386 //----------------------------------------------------------------------------
01387 int fei::VectorSpace::getGlobalNumBlkIndices() const
01388 {
01389   int numBlkIndices = 0;
01390   unsigned len = globalIDOffsets_.size();
01391   if (len > 0) {
01392     numBlkIndices = globalIDOffsets_[len-1];
01393   }
01394 
01395   return(numBlkIndices);
01396 }
01397 
01398 //----------------------------------------------------------------------------
01399 int fei::VectorSpace::getBlkIndices_Owned(int lenBlkIndices,
01400                                           int* globalBlkIndices,
01401                                           int* blkSizes,
01402                                           int& numBlkIndices)
01403 {
01404   if (!sharedRecordsSynchronized_) {
01405     numBlkIndices = 0;
01406     return(-1);
01407   }
01408 
01409   int localProc = fei::localProc(comm_);
01410   fei::BlkIndexAccessor blkIndAccessor(localProc, lenBlkIndices,
01411                                            globalBlkIndices, blkSizes);
01412   runRecords(blkIndAccessor);
01413 
01414   numBlkIndices = blkIndAccessor.numBlkIndices_;
01415 
01416   return(0);
01417 }
01418 
01419 //----------------------------------------------------------------------------
01420 int fei::VectorSpace::getRecordCollection(int idType,
01421                                           snl_fei::RecordCollection*& records)
01422 {
01423   int idx = fei::binarySearch(idType, idTypes_);
01424   if (idx < 0) return(-1);
01425 
01426   records = recordCollections_[idx];
01427   return(0);
01428 }
01429 
01430 //----------------------------------------------------------------------------
01431 int fei::VectorSpace::getRecordCollection(int idType,
01432                                           const snl_fei::RecordCollection*& records) const
01433 {
01434   int idx = fei::binarySearch(idType, idTypes_);
01435   if (idx < 0) return(-1);
01436 
01437   records = recordCollections_[idx];
01438   return(0);
01439 }
01440 
01441 //----------------------------------------------------------------------------
01442 void fei::VectorSpace::setOwners_lowestSharing()
01443 {
01444   //first, add localProc to each of the sharing-proc lists, in case it wasn't
01445   //included via initSharedIDs().
01446   int localProc = fei::localProc(comm_);
01447 
01448   std::map<int, fei::SharedIDs<int> >::iterator
01449     s_iter = sharedIDTables_.begin(), s_end = sharedIDTables_.end();
01450 
01451   for(; s_iter != s_end; ++s_iter) {
01452     fei::SharedIDs<int>::map_type& shid_table = s_iter->second.getSharedIDs();
01453 
01454     fei::SharedIDs<int>::map_type::iterator
01455       t_iter = shid_table.begin(),
01456       t_end = shid_table.end();
01457 
01458     for(; t_iter != t_end; ++t_iter) {
01459       std::set<int>& shProcs = t_iter->second;
01460       shProcs.insert(localProc);
01461     }
01462   }
01463 
01464   //now set the owningProcs for the SharedIDs records, and the owning procs on
01465   //the appropriate records in the recordCollections. Set the owner to be the
01466   //lowest-numbered sharing proc in all cases.
01467   s_iter = sharedIDTables_.begin();
01468   for(; s_iter != s_end; ++s_iter) {
01469     fei::SharedIDs<int>::map_type& shid_table = s_iter->second.getSharedIDs();
01470 
01471     std::vector<int>& owningProcs = s_iter->second.getOwningProcs();
01472 
01473     int len = shid_table.size();
01474     owningProcs.resize(len);
01475     int j = 0;
01476 
01477     fei::SharedIDs<int>::map_type::iterator
01478       t_iter = shid_table.begin(),
01479       t_end = shid_table.end();
01480 
01481     for(; t_iter != t_end; ++t_iter, ++j) {
01482       std::set<int>& shProcs = (*t_iter).second;
01483       int lowest = *(shProcs.begin());
01484       owningProcs[j] = lowest;
01485     }
01486     
01487     int idx = fei::binarySearch(s_iter->first, idTypes_);
01488     if (idx < 0) {
01489       throw std::runtime_error("internal error in fei::VectorSapce::setOwners_lowestSharing");
01490     }
01491 
01492     if (output_level_ >= fei::FULL_LOGS && output_stream_ != NULL) {
01493       recordCollections_[idx]->setDebugOutput(output_stream_);
01494     }
01495 
01496     recordCollections_[idx]->setOwners_lowestSharing(s_iter->second);
01497   }
01498 }
01499 
01500 //----------------------------------------------------------------------------
01501 int fei::VectorSpace::calculateGlobalIndices()
01502 {
01503   int localProc = fei::localProc(comm_);
01504   int numProcs = fei::numProcs(comm_);
01505   std::vector<int> localOffsets(numProcs +1, 0);
01506   std::vector<int> localIDOffsets(numProcs+1, 0);
01507   globalOffsets_.resize(numProcs + 1, 0);
01508 
01509   globalIDOffsets_.assign(globalIDOffsets_.size(), 0);
01510 
01511   //first we'll calculate the number of local degrees of freedom, and the
01512   //number of local identifiers.
01513   fei::RecordAttributeCounter counter(localProc);
01514   runRecords(counter);
01515 
01516   int numLocalDOF = counter.numLocalDOF_;
01517   int numLocalIDs = counter.numLocalIDs_;
01518   int numRemoteSharedDOF = counter.numRemoteSharedDOF_;
01519 
01520   eqnNumbers_.resize(numLocalDOF+numRemoteSharedDOF);
01521 
01522   localOffsets[localProc] = numLocalDOF;
01523   CHK_MPI( fei::GlobalMax(comm_, localOffsets, globalOffsets_) );
01524 
01525   localIDOffsets[localProc] = numLocalIDs;
01526   CHK_MPI( fei::GlobalMax(comm_, localIDOffsets, globalIDOffsets_) );
01527 
01528   //Now globalOffsets_ contains numLocalDOF for proc i, in the i-th position.
01529   //(and similarly for globalIDOffsets_)
01530   //So all we need to do is turn that data into global-offsets (i.e., the
01531   //starting global-offset for each processor).
01532   int localOffset = 0;
01533   int localIDOffset = 0;
01534   for(int p=0; p<numProcs; ++p) {
01535     numLocalDOF = globalOffsets_[p];
01536     globalOffsets_[p] = localOffset;
01537     localOffset += numLocalDOF;
01538     numLocalIDs = globalIDOffsets_[p];
01539     globalIDOffsets_[p] = localIDOffset;
01540     localIDOffset += numLocalIDs;
01541   }
01542   globalOffsets_[numProcs] = localOffset;
01543   globalIDOffsets_[numProcs] = localIDOffset;
01544 
01545   firstLocalOffset_ = globalOffsets_[localProc];
01546   lastLocalOffset_ = globalOffsets_[localProc+1] - 1;
01547 
01548   if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
01549     FEI_OSTREAM& os = *output_stream_;
01550     os <<dbgprefix_<< "  firstLocalOffset_: " << firstLocalOffset_ << ", lastLocalOffset_: "
01551       << lastLocalOffset_ << FEI_ENDL;
01552   }
01553 
01554   //Now we're ready to set the equation-numbers on all local records.
01555   CHK_ERR( setLocalEqnNumbers() );
01556 
01557   return(0);
01558 }
01559 
01560 //----------------------------------------------------------------------------
01561 int fei::VectorSpace::synchronizeSharedRecords()
01562 {
01563   if (output_level_ >= fei::FULL_LOGS && output_stream_ != NULL) {
01564     FEI_OSTREAM& os = *output_stream_;
01565     os<<dbgprefix_<<"#synchronizeSharedRecords num-field-masks: "<<fieldMasks_.size()<<FEI_ENDL;
01566     for(unsigned fm=0; fm<fieldMasks_.size(); ++fm) {
01567       os << dbgprefix_<<"#     maskID["<<fm<<"]: " << fieldMasks_[fm]->getMaskID() << FEI_ENDL;
01568     }
01569   }
01570 
01571   bool safetyCheck = checkSharedIDs_;
01572 
01573   int numProcs = fei::numProcs(comm_);
01574   int localProc = fei::localProc(comm_);
01575   int numShTables = sharedIDTables_.size();
01576 
01577   if (numProcs < 2) {
01578     sharedRecordsSynchronized_ = true;
01579     return(0);
01580   }
01581 
01582   if (safetyCheck) {
01583     if (localProc == 0 && output_level_ >= fei::BRIEF_LOGS) {
01584       FEI_COUT << "fei::VectorSpace: global consistency-check of shared ID"
01585            << " data (involves all-to-all communication). This is done "
01586            << "only if requested by parameter 'FEI_CHECK_SHARED_IDS true'."
01587            << FEI_ENDL;
01588     }
01589 
01590     int totalNumShTables = 0;
01591     CHK_ERR( fei::GlobalSum(comm_, numShTables, totalNumShTables) );
01592     if (totalNumShTables != numShTables * numProcs) {
01593       //not all processors have the same number of shared-id tables. This means
01594       //that one or more processors is 'empty', which may not be an error. But
01595       //it does mean that the safety check can't be performed because the safety
01596       //check involves all-to-all communication and if one or more processors
01597       //don't participate then we'll hang...
01598       safetyCheck = false;
01599     }
01600   }
01601 
01602   //Create a list of fei::comm_maps which will be the communication-patterns
01603   //for each of the sharedIDTables. A sharedIDTable maps lists of processors
01604   //to each shared ID. The communication-pattern will be a transpose of
01605   //that, mapping lists of IDs to owning or sharing processors.
01606   //
01607 
01608   int local_err = 0;
01609 
01610   for(size_t i=0; i<idTypes_.size(); ++i) {
01611 
01612     fei::SharedIDs<int>& shared = getSharedIDs_private(idTypes_[i]);
01613     
01614     //now create/initialize ownerPatterns which map owning processors to lists
01615     //of ids that are shared locally, and sharerPatterns which map sharing
01616     //processors to lists of ids that are owned locally.
01617     fei::comm_map* ownerPattern = new fei::comm_map(1, numProcs);
01618 
01619     ownerPatterns_[idTypes_[i]] = ownerPattern;
01620 
01621     fei::comm_map* sharerPattern = 
01622             new fei::comm_map(1, numProcs);
01623 
01624     sharerPatterns_[idTypes_[i]] = sharerPattern;
01625 
01626     std::vector<int>& owningProcs = shared.getOwningProcs();
01627 
01628     fei::SharedIDs<int>::map_type& shtable = shared.getSharedIDs();
01629 
01630     fei::SharedIDs<int>::map_type::iterator
01631       s_iter = shtable.begin(),
01632       s_end = shtable.end();
01633     
01634     int j = 0;
01635     for(; s_iter != s_end; ++s_iter, ++j) {
01636       int ID = (*s_iter).first;
01637       std::set<int>& shProcs = (*s_iter).second;
01638 
01639       int owner = owningProcs[j];
01640 
01641       if (owner == localProc) {
01642         std::set<int>::const_iterator
01643           p_iter = shProcs.begin(),
01644           p_end = shProcs.end();
01645         for(; p_iter != p_end; ++p_iter) {
01646           if (*p_iter != localProc) {
01647             sharerPattern->addIndices(*p_iter, 1, &ID);
01648           }
01649         }
01650       }
01651       else {
01652         ownerPattern->addIndices(owner, 1, &ID);
01653       }
01654     }
01655 
01656     if (safetyCheck) {
01657       fei::comm_map* checkPattern = NULL;
01658       CHK_ERR( fei::mirrorCommPattern(comm_, ownerPattern, checkPattern) );
01659       fei::Barrier(comm_);
01660       if (output_level_ >= fei::FULL_LOGS && output_stream_ != NULL) {
01661         FEI_OSTREAM& os = *output_stream_;
01662         fei::comm_map::map_type& owner_map = ownerPattern->getMap();
01663         int numKeys = owner_map.size();
01664         fei::comm_map::map_type::const_iterator
01665           omap_iter = owner_map.begin(),
01666           omap_end = owner_map.end();
01667 
01668         os << dbgprefix_<<"#  synchronizeSharedRecords" << FEI_ENDL
01669            << dbgprefix_<<"#  ownerPattern, num-procs-to-send-to: " << numKeys << FEI_ENDL;
01670         for(int sk=0; omap_iter != omap_end; ++sk, ++omap_iter) {
01671           os << dbgprefix_<<"#    sendProc["<<sk<<"]: " << omap_iter->first << " IDs: ";
01672           fei::comm_map::row_type::const_iterator
01673             val_iter = omap_iter->second->begin(),
01674             val_end =  omap_iter->second->end();
01675           for(; val_iter != val_end; ++val_iter) {
01676             os << *val_iter << " ";
01677           }
01678           os << FEI_ENDL;
01679         }
01680 
01681         fei::comm_map::map_type& check_map = checkPattern->getMap();
01682         int numCheckKeys = check_map.size();
01683         fei::comm_map::map_type::const_iterator
01684           cmap_iter = check_map.begin(),
01685           cmap_end = check_map.end();
01686 
01687         os <<dbgprefix_<< "#  synchronizeSharedRecords" << FEI_ENDL
01688            <<dbgprefix_<< "#  checkPattern (send mirror), num-procs: "
01689            << numCheckKeys << FEI_ENDL;
01690         for(int sk=0; cmap_iter != cmap_end; ++sk, ++cmap_iter) {
01691           os <<dbgprefix_<< "#    proc["<<sk<<"]: " << cmap_iter->first << " IDs: ";
01692           fei::comm_map::row_type::const_iterator
01693           val_iter = cmap_iter->second->begin(),
01694             val_end = cmap_iter->second->end();
01695           for(; val_iter != val_end; ++val_iter) {
01696             os << *val_iter << " ";
01697           }
01698           os << FEI_ENDL;
01699         }
01700       }
01701 
01702       int err = 0;
01703       bool quiet = false;
01704       if (!checkPattern->equal(*sharerPattern, quiet)) {
01705         //FEI_CERR << "shared-ID safety-check failed..."
01706         //     << FEI_ENDL;
01707         err = -1;
01708       }
01709       int globalErr = 0;
01710       CHK_ERR( fei::GlobalSum(comm_, err, globalErr) );
01711 
01712       delete checkPattern;
01713 
01714       if (globalErr != 0) {
01715         return(globalErr);
01716       }
01717     }
01718 
01719     local_err += exchangeFieldInfo(ownerPattern, sharerPattern,
01720                                    recordCollections_[i], fieldMasks_);
01721   }
01722 
01723   int global_err = 0;
01724   CHK_ERR( fei::GlobalSum(comm_, local_err, global_err) );
01725 
01726   if (global_err != 0) {
01727     ERReturn(-1);
01728   }
01729 
01730   sharedRecordsSynchronized_ = true;
01731 
01732   return(0);
01733 }
01734 
01735 //----------------------------------------------------------------------------
01736 int fei::VectorSpace::exchangeGlobalIndices()
01737 {
01738   for(size_t i=0; i<idTypes_.size(); ++i) {
01739     snl_fei::RecordMsgHandler recmsghndlr(fei::localProc(comm_),
01740                                  recordCollections_[i], *ptBlkMap_,
01741                                           fieldMasks_, eqnNumbers_);
01742     recmsghndlr.setTask(snl_fei::RecordMsgHandler::_EqnNumbers_);
01743 
01744     recmsghndlr.setSendPattern(sharerPatterns_[idTypes_[i]]);
01745     recmsghndlr.setRecvPattern(ownerPatterns_[idTypes_[i]]);
01746     CHK_ERR( fei::exchange(comm_, &recmsghndlr) );
01747   }
01748 
01749   return(0);
01750 }
01751 
01752 //----------------------------------------------------------------------------
01753 void fei::VectorSpace::runRecords(fei::Record_Operator<int>& record_op)
01754 {
01755   for(size_t rec=0; rec<recordCollections_.size(); ++rec) {
01756     snl_fei::RecordCollection* records = recordCollections_[rec];
01757 
01758     snl_fei::RecordCollection::map_type& rmap = records->getRecords();
01759 
01760     snl_fei::RecordCollection::map_type::iterator
01761       r_iter = rmap.begin(),
01762       r_end  = rmap.end();
01763 
01764     for(; r_iter != r_end; ++r_iter) {
01765       fei::Record<int>& thisrecord = *((*r_iter).second);
01766 
01767       record_op(thisrecord);
01768     }
01769   }
01770 }
01771 
01772 //----------------------------------------------------------------------------
01773 int fei::VectorSpace::setLocalEqnNumbers()
01774 {
01775   int proc = fei::localProc(comm_);
01776   int eqnNumber = firstLocalOffset_;
01777   int idNumber = globalIDOffsets_[proc];
01778 
01779   int numProcs = fei::numProcs(comm_);
01780   int localProc = fei::localProc(comm_);
01781 
01782   if (ptBlkMap_ != NULL) {
01783     delete ptBlkMap_;
01784   }
01785   ptBlkMap_ = new snl_fei::PointBlockMap;
01786 
01787   int maxNumIndices = 0;
01788   for(unsigned i=0; i<fieldMasks_.size(); ++i) {
01789     if (fieldMasks_[i]->getNumIndices() > maxNumIndices) {
01790       maxNumIndices = fieldMasks_[i]->getNumIndices();
01791     }
01792   }
01793 
01794   if (maxNumIndices == 1) {
01795     ptBlkMap_->setPtEqualBlk();
01796   }
01797 
01798   FEI_OSTREAM* id2eqnStream = NULL;
01799   if (output_level_ >= fei::BRIEF_LOGS) {
01800     std::string path = fei::LogManager::getLogManager().getOutputPath();
01801     if (path == "") path = ".";
01802     FEI_OSTRINGSTREAM osstr;
01803     osstr << path << "/fei_ID2Eqn";
01804     if (name_.size() > 0) osstr << "_" << name_;
01805     osstr << "." <<numProcs<< "." << localProc;
01806 
01807     id2eqnStream = new FEI_OFSTREAM(osstr.str().c_str(), IOS_OUT);
01808     FEI_OSTREAM& os = *id2eqnStream;
01809     os << "# Each line contains:\n#   ID   blk-eqn   eqn" << FEI_ENDL;
01810   }
01811 
01812   int eqnNumberOffset = 0;
01813   int maxNumDOF = 0;
01814 
01815   for(size_t rec=0; rec<recordCollections_.size(); ++rec) {
01816     snl_fei::RecordCollection* records = recordCollections_[rec];
01817 
01818     snl_fei::RecordCollection::map_type& rmap = records->getRecords();
01819 
01820     snl_fei::RecordCollection::map_type::iterator
01821       r_iter = rmap.begin(),
01822       r_end = rmap.end();
01823 
01824     int* eqnNumPtr = &eqnNumbers_[0];
01825 
01826     for(; r_iter != r_end; ++r_iter) {
01827       fei::Record<int>& thisrecord = *((*r_iter).second);
01828 
01829       fei::FieldMask* mask = thisrecord.getFieldMask();
01830 
01831       thisrecord.setOffsetIntoEqnNumbers(eqnNumberOffset);
01832 
01833       int owner = thisrecord.getOwnerProc();
01834       if (owner == proc) {
01835         thisrecord.setNumber(idNumber++);
01836       }
01837 
01838       int* eqnNumbers = eqnNumPtr
01839                       + thisrecord.getOffsetIntoEqnNumbers();
01840 
01841       int numDOF = mask->getNumIndices();
01842       eqnNumberOffset += numDOF;
01843 
01844       if (output_level_ >= fei::BRIEF_LOGS) {
01845         for(int ii=0; ii<numDOF; ++ii) {
01846           if (isLogEqn(eqnNumber+ii) && output_stream_ != NULL) {
01847             FEI_OSTREAM& os = *output_stream_;
01848             os <<dbgprefix_<< "setLocalEqnNumbers: ID "<<thisrecord.getID()
01849                << " <--> eqn " << eqnNumber+ii<<FEI_ENDL;
01850           }
01851         }
01852       }
01853 
01854       if (owner == proc) {
01855         int offset = 0;
01856         for(int n=0; n<numDOF; ++n) {
01857           eqnNumbers[offset++] = eqnNumber++;
01858         }
01859       }
01860 
01861       if (numDOF > maxNumDOF) maxNumDOF = numDOF;
01862 
01863       if (owner == proc) {
01864         int thiseqn = eqnNumber-numDOF;
01865         int thisrecordnumber = thisrecord.getNumber();
01866         if (maxNumIndices > 1) {
01867           CHK_ERR( ptBlkMap_->setEqn(thiseqn, thisrecordnumber, numDOF) );
01868           if (numDOF > 1) {
01869             for(int i=1; i<numDOF; ++i) {
01870               CHK_ERR( ptBlkMap_->setEqn(thiseqn+i, thisrecordnumber, numDOF) );
01871             }
01872           }
01873         }
01874       }
01875 
01876       if (id2eqnStream != NULL) {
01877         if (owner == proc) {
01878           FEI_OSTREAM& os = *id2eqnStream;
01879           for(int n=0; n<numDOF; ++n) {
01880             os << thisrecord.getID() << " " << thisrecord.getNumber() << " "
01881                << eqnNumber-numDOF+n<<FEI_ENDL;
01882           }
01883         }
01884       }
01885     }
01886 
01887   }
01888 
01889   ptBlkMap_->setMaxBlkEqnSize(maxNumDOF);
01890 
01891   int globalMaxNumDOF;
01892   CHK_ERR( fei::GlobalMax(comm_, maxNumDOF, globalMaxNumDOF) );
01893 
01894   if (globalMaxNumDOF == 1) {
01895     ptBlkMap_->setPtEqualBlk();
01896   }
01897 
01898   delete id2eqnStream;
01899 
01900   return(0);
01901 }
01902 
01903 //----------------------------------------------------------------------------
01904 int fei::VectorSpace::exchangeFieldInfo(fei::comm_map* ownerPattern,
01905                                         fei::comm_map* sharerPattern,
01906                                   snl_fei::RecordCollection* recordCollection,
01907                                       std::vector<fei::FieldMask*>& fieldMasks)
01908 {
01909   //ownerPattern maps owning processors to lists of IDs that we (the local
01910   //processor) share.
01911   //
01912   //sharerPattern maps sharing processors to lists of IDs that we own.
01913   //
01914   //In this function we need to perform these tasks:
01915   //1. processors exchange and combine their sets of field-masks so that all
01916   //   processors will have the super-set of field-masks that they need.
01917   //2. sharing processors send maskIDs for their shared IDs to the owning
01918   //   processors. The owning processors then combine masks as necessary to
01919   //   make sure each shared ID has the union of field-masks that are held by
01920   //   each of the sharing processors. all owning processors now send their
01921   //   maskIDs out to the sharing processors. At this point all processors
01922   //   should have the same view of the masks that belong on shared IDs.
01923   //3. exchange info describing inactive DOF offsets for shared records.
01924   //
01925 
01926   int numProcs = fei::numProcs(comm_);
01927   if (numProcs < 2) return(0);
01928 
01929   if (ptBlkMap_ == 0) ptBlkMap_ = new snl_fei::PointBlockMap;
01930 
01931   snl_fei::RecordMsgHandler recMsgHandler(fei::localProc(comm_),
01932                                           recordCollection,
01933                                           *ptBlkMap_,
01934                                           fieldMasks, eqnNumbers_);
01935 
01936   //Step 1a.
01937   recMsgHandler.setTask(snl_fei::RecordMsgHandler::_FieldMasks_);
01938 
01939   recMsgHandler.setSendPattern(ownerPattern);
01940   recMsgHandler.setRecvPattern(sharerPattern);
01941   CHK_ERR( fei::exchange(comm_, &recMsgHandler) );
01942 
01943   //Step 2a.
01944   recMsgHandler.setTask(snl_fei::RecordMsgHandler::_MaskIDs_);
01945 
01946   recMsgHandler.setSendPattern(ownerPattern);
01947   recMsgHandler.setRecvPattern(sharerPattern);
01948   CHK_ERR( fei::exchange(comm_, &recMsgHandler) );
01949 
01950   //Step 1b.
01951   recMsgHandler.setTask(snl_fei::RecordMsgHandler::_FieldMasks_);
01952 
01953   recMsgHandler.setSendPattern(sharerPattern);
01954   recMsgHandler.setRecvPattern(ownerPattern);
01955   CHK_ERR( fei::exchange(comm_, &recMsgHandler) );
01956 
01957   //Step 2b.
01958   recMsgHandler.setTask(snl_fei::RecordMsgHandler::_MaskIDs_);
01959 
01960   recMsgHandler.setSendPattern(sharerPattern);
01961   recMsgHandler.setRecvPattern(ownerPattern);
01962   CHK_ERR( fei::exchange(comm_, &recMsgHandler) );
01963 
01964   return(0);
01965 }
01966 
01967 //----------------------------------------------------------------------------
01968 fei::FieldDofMap<int>&
01969 fei::VectorSpace::getFieldDofMap()
01970 {
01971   return fieldDofMap_;
01972 }
01973 
01974 //----------------------------------------------------------------------------
01975 void fei::VectorSpace::setName(const char* name)
01976 {
01977   if (name == NULL) return;
01978 
01979   if (name_ == name) return;
01980 
01981   name_ = name;
01982   dbgprefix_ = "VecSpc_"+name_+": ";
01983 }
01984 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends
Generated on Wed Apr 13 10:08:24 2011 for FEI by  doxygen 1.6.3