FEI Version of the Day
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   int localNeedToCompute = need_to_compute_shared_ids ? 1 : 0;
00628   int globalNeedToCompute = 0;
00629   CHK_ERR( fei::GlobalMin(comm_, localNeedToCompute, globalNeedToCompute) );
00630   need_to_compute_shared_ids = globalNeedToCompute==1 ? true : false;
00631   if (need_to_compute_shared_ids) {
00632     compute_shared_ids();
00633   }
00634 
00635   //setOwners_lowestSharing is a local operation (no communication), which
00636   //assumes that each processor holds CORRECT (globally symmetric)
00637   //shared-id/sharing-proc tables. No correctness-checking is performed here.
00638   setOwners_lowestSharing();
00639 
00640   //synchronizeSharedRecords ensures that each sharing processor has the same
00641   //view of the shared records, with respect to the layout of fields, which
00642   //determines how many DOFs and equation-numbers reside at each ID.
00643   //This involves inter-processor communication.
00644   if ( synchronizeSharedRecords() != 0) {
00645     return(-1);
00646   }
00647 
00648   //calculateGlobalIndices is also a global operation.
00649   if ( calculateGlobalIndices() != 0) {
00650     return(-1);
00651   }
00652 
00653   //finally we need to exchange global indices for shared records. i.e., the
00654   //processors that own shared records, need to send the global indices for
00655   //those records to the sharing-but-not-owning processors.
00656   if (fei::numProcs(comm_) > 1) {
00657     if ( exchangeGlobalIndices() != 0) {
00658       return(-1);
00659     }
00660   }
00661 
00662   newInitData_ = false;
00663 
00664   return(0);
00665 }
00666 
00667 //----------------------------------------------------------------------------
00668 int fei::VectorSpace::getGlobalIndex(int idType,
00669                                      int ID,
00670                                      int fieldID,
00671                                      int fieldOffset,
00672                                      int whichComponentOfField,
00673                                      int& globalIndex)
00674 {
00675   int idindex = fei::binarySearch(idType, idTypes_);
00676   if (idindex < 0) return(-1);
00677 
00678   unsigned fieldSize = 0;
00679   if (fieldOffset > 0) {
00680     fieldSize = getFieldSize(fieldID);
00681   }
00682 
00683   try {
00684     globalIndex = recordCollections_[idindex]->getGlobalIndex(ID,
00685                                                          fieldID,
00686                                                          fieldSize,
00687                                                          fieldOffset,
00688                                                          whichComponentOfField,
00689                                                          &eqnNumbers_[0]);
00690   }
00691   catch (std::runtime_error& exc) {
00692     FEI_OSTRINGSTREAM osstr;
00693     osstr << "VectorSpace::getGlobalIndex caught exception: " << exc.what();
00694     FEI_CERR << osstr.str()<<FEI_ENDL;
00695     ERReturn(-1);
00696   }
00697 
00698   return(0);
00699 }
00700 
00701 //----------------------------------------------------------------------------
00702 int fei::VectorSpace::getGlobalIndex(int idType,
00703                                      int ID,
00704                                      int fieldID,
00705                                      int& globalIndex)
00706 {
00707   return( getGlobalIndex(idType, ID, fieldID, 0, 0, globalIndex) );
00708 }
00709 
00710 //----------------------------------------------------------------------------
00711 int fei::VectorSpace::getGlobalBlkIndex(int idType,
00712                                         int ID,
00713                                         int& globalBlkIndex)
00714 {
00715   int idindex = fei::binarySearch(idType, idTypes_);
00716   if (idindex < 0) return(-1);
00717 
00718   CHK_ERR(recordCollections_[idindex]->getGlobalBlkIndex(ID, globalBlkIndex));
00719 
00720   return(0);
00721 }
00722 
00723 //----------------------------------------------------------------------------
00724 int fei::VectorSpace::getGlobalIndices(int numIDs,
00725                                        const int* IDs,
00726                                        int idType,
00727                                        int fieldID,
00728                                        int* globalIndices)
00729 {
00730   int idindex = fei::binarySearch(idType, idTypes_);
00731   if (idindex < 0) return(-1);
00732 
00733   unsigned fieldSize = getFieldSize(fieldID);
00734   int offset = 0;
00735 
00736   for(int i=0; i<numIDs; ++i) {
00737     try {
00738       globalIndices[offset] =
00739         recordCollections_[idindex]->getGlobalIndex(IDs[i],
00740                                                     fieldID,
00741                                                     fieldSize,
00742                                                     0, 0,
00743                                                     &eqnNumbers_[0]);
00744       if (fieldSize>1) {
00745         int eqn = globalIndices[offset];
00746         for(unsigned j=1; j<fieldSize; ++j) {
00747           globalIndices[offset+j] = eqn+j;
00748         }
00749       }
00750     }
00751     catch (...) {
00752       for(unsigned j=0; j<fieldSize; ++j) {
00753         globalIndices[offset+j] = -1;
00754       }
00755     }
00756 
00757     offset += fieldSize;
00758   }
00759 
00760   return(0);
00761 }
00762 
00763 //----------------------------------------------------------------------------
00764 int fei::VectorSpace::getGlobalBlkIndices(int numIDs,
00765                                           const int* IDs,
00766                                           int idType,
00767                                           int* globalBlkIndices)
00768 {
00769   int err;
00770   int idindex = fei::binarySearch(idType, idTypes_);
00771   if (idindex < 0) return(-1);
00772 
00773   int offset = 0;
00774 
00775   for(int i=0; i<numIDs; ++i) {
00776     err = recordCollections_[idindex]->getGlobalBlkIndex(IDs[i],
00777                                                      globalBlkIndices[offset]);
00778     if (err != 0) {
00779       globalBlkIndices[offset] = -1;
00780     }
00781     ++offset;
00782   }
00783 
00784   return(0);
00785 }
00786 
00787 //----------------------------------------------------------------------------
00788 int fei::VectorSpace::getGlobalIndices(int numIDs,
00789                                        const int* IDs,
00790                                        const int* idTypes,
00791                                        const int* fieldIDs,
00792                                        int* globalIndices)
00793 {
00794   int err;
00795   int offset = 0;
00796   for(int i=0; i<numIDs; ++i) {
00797     unsigned fieldSize = getFieldSize(fieldIDs[i]);
00798     err = getGlobalIndex(idTypes[i], IDs[i], fieldIDs[i], 0, 0,
00799                          globalIndices[offset]);
00800     if (err) {
00801       for(unsigned j=1; j<fieldSize; ++j) {
00802         globalIndices[offset+j] = -1;
00803       }
00804     }
00805     else {
00806       if (fieldSize>1) {
00807         int eqn = globalIndices[offset];
00808         for(unsigned j=1; j<fieldSize; ++j) {
00809           globalIndices[offset+j] = eqn+j;
00810         }
00811       }
00812     }
00813     offset += fieldSize;
00814   }
00815 
00816   return(0);
00817 }
00818 
00819 //----------------------------------------------------------------------------
00820 void fei::VectorSpace::getGlobalBlkIndices(const fei::Pattern* pattern,
00821                                            const fei::Record<int>*const* records,
00822                                            std::vector<int>& indices)
00823 {
00824   int numRecords = pattern->getNumIDs();
00825   indices.resize(numRecords);
00826   int numIndices;
00827   getGlobalBlkIndices(numRecords, records, numRecords, &indices[0],
00828                       numIndices);
00829 }
00830 
00831 //----------------------------------------------------------------------------
00832 void fei::VectorSpace::getGlobalIndices(const fei::Pattern* pattern,
00833                                         const fei::Record<int>*const* records,
00834                                         std::vector<int>& indices)
00835 {
00836   int numRecords = pattern->getNumIDs();
00837   int numIndices = pattern->getNumIndices();
00838   indices.resize(numIndices);
00839   int* indices_ptr = &indices[0];
00840 
00841   fei::Pattern::PatternType pType = pattern->getPatternType();
00842 
00843   if (pType == fei::Pattern::GENERAL ||
00844       pType == fei::Pattern::SINGLE_IDTYPE) {
00845     const int* numFieldsPerID = pattern->getNumFieldsPerID();
00846     const int* fieldIDs = pattern->getFieldIDs();
00847     int totalNumFields = pattern->getTotalNumFields();
00848 
00849     std::vector<int> fieldSizes(totalNumFields);
00850 
00851     for(int j=0; j<totalNumFields; ++j) {
00852       fieldSizes[j] = getFieldSize(fieldIDs[j]);
00853     }
00854 
00855     getGlobalIndices(numRecords, records, numFieldsPerID,
00856                      fieldIDs, &(fieldSizes[0]),
00857                      numIndices, indices_ptr, numIndices);
00858   }
00859   else if (pType == fei::Pattern::SIMPLE) {
00860     const int* fieldIDs = pattern->getFieldIDs();
00861 
00862     int fieldID = fieldIDs[0];
00863     unsigned fieldSize = getFieldSize(fieldID);
00864 
00865     getGlobalIndices(numRecords, records,
00866                      fieldID, fieldSize,
00867                      numIndices, indices_ptr, numIndices);
00868   }
00869   else if (pType == fei::Pattern::NO_FIELD) {
00870     getGlobalBlkIndices(numRecords, records, numIndices, indices_ptr, numIndices);
00871   }
00872 }
00873 
00874 //----------------------------------------------------------------------------
00875 void fei::VectorSpace::getGlobalIndices(int numRecords,
00876                                         const fei::Record<int>*const* records,
00877                                         int fieldID,
00878                                         int fieldSize,
00879                                         int indicesAllocLen,
00880                                         int* indices,
00881                                         int& numIndices)
00882 {
00883   numIndices = 0;
00884   int* eqnPtr = &eqnNumbers_[0];
00885 
00886   int len = numRecords;
00887   if (len*fieldSize >= indicesAllocLen) {
00888     len = indicesAllocLen/fieldSize;
00889   }
00890 
00891   if (fieldSize == 1 && simpleProblem_) {
00892     for(int i=0; i<len; ++i) {
00893       const fei::Record<int>* record = records[i];
00894       indices[numIndices++] = *(eqnPtr+record->getOffsetIntoEqnNumbers());
00895     }    
00896     return;
00897   }
00898 
00899   if (fieldSize == 1) {
00900     for(int i=0; i<len; ++i) {
00901       const fei::Record<int>* record = records[i];
00902 
00903       int eqnOffset = 0, numInstances = 0;
00904       int* eqnNumbers = eqnPtr+record->getOffsetIntoEqnNumbers();
00905 
00906       const fei::FieldMask* fieldMask = record->getFieldMask();
00907       fieldMask->getFieldEqnOffset(fieldID, eqnOffset, numInstances);
00908       indices[numIndices++] = eqnNumbers[eqnOffset];
00909     }
00910   }
00911   else {
00912     for(int i=0; i<len; ++i) {
00913       const fei::Record<int>* record = records[i];
00914 
00915       int* eqnNumbers = eqnPtr+record->getOffsetIntoEqnNumbers();
00916 
00917       int eqnOffset = 0, numInstances = 0;
00918       if (!simpleProblem_) {
00919         const fei::FieldMask* fieldMask = record->getFieldMask();
00920         fieldMask->getFieldEqnOffset(fieldID, eqnOffset, numInstances);
00921       }
00922 
00923       for(int fs=0; fs<fieldSize; ++fs) {
00924         indices[numIndices++] = eqnNumbers[eqnOffset+fs];
00925       }
00926     }
00927   }
00928 }
00929 
00930 //----------------------------------------------------------------------------
00931 void fei::VectorSpace::getGlobalIndices(int numRecords,
00932                                         const fei::Record<int>*const* records,
00933                                         const int* numFieldsPerID,
00934                                         const int* fieldIDs,
00935                                         const int* fieldSizes,
00936                                         int indicesAllocLen,
00937                                         int* indices,
00938                                         int& numIndices)
00939 {
00940   numIndices = 0;
00941   int fld_offset = 0;
00942   int* eqnPtr = &eqnNumbers_[0];
00943 
00944   for(int i=0; i<numRecords; ++i) {
00945     const fei::Record<int>* record = records[i];
00946 
00947     const fei::FieldMask* fieldMask = record->getFieldMask();
00948     int* eqnNumbers = eqnPtr + record->getOffsetIntoEqnNumbers();
00949 
00950     for(int nf=0; nf<numFieldsPerID[i]; ++nf) {
00951       int eqnOffset = 0, numInstances = 0;
00952       if (!simpleProblem_) {
00953         fieldMask->getFieldEqnOffset(fieldIDs[fld_offset], eqnOffset, numInstances);
00954       }
00955 
00956       for(int fs=0; fs<fieldSizes[fld_offset]; ++fs) {
00957         indices[numIndices++] = eqnNumbers[eqnOffset+fs];
00958       }
00959 
00960       ++fld_offset;
00961     }
00962   }
00963 }
00964 
00965 //----------------------------------------------------------------------------
00966 void fei::VectorSpace::getGlobalBlkIndices(int numRecords,
00967                                            const fei::Record<int>*const* records,
00968                                            int indicesAllocLen,
00969                                            int* indices,
00970                                            int& numIndices)
00971 {
00972   numIndices = 0;
00973   for(int i=0; i<numRecords; ++i) {
00974     if (numIndices < indicesAllocLen) {
00975       indices[numIndices++] = records[i]->getNumber();
00976     }
00977     else ++numIndices;
00978   }
00979 }
00980 
00981 //----------------------------------------------------------------------------
00982 int fei::VectorSpace::getGlobalIndex(int idType,
00983                                      int ID,
00984                                      int& globalIndex)
00985 {
00986   int idindex = fei::binarySearch(idType, idTypes_);
00987   if (idindex < 0) return(-1);
00988 
00989   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
00990   if (record == NULL) {
00991     ERReturn(-1);
00992   }
00993 
00994   const int* eqnNums = &eqnNumbers_[0]
00995                      + record->getOffsetIntoEqnNumbers();
00996 
00997   if (eqnNums != NULL) { globalIndex = eqnNums[0]; return(0); }
00998   return(-1);
00999 }
01000 
01001 //----------------------------------------------------------------------------
01002 int fei::VectorSpace::getNumDegreesOfFreedom(int idType,
01003                                              int ID)
01004 {
01005   int idindex = fei::binarySearch(idType, idTypes_);
01006   if (idindex < 0) return(0);
01007 
01008   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
01009   if (record == NULL) {
01010     ERReturn(-1);
01011   }
01012 
01013   return( record->getFieldMask()->getNumIndices() );
01014 }
01015 
01016 //----------------------------------------------------------------------------
01017 int fei::VectorSpace::getNumFields()
01018 {
01019   return(fieldDatabase_.size());
01020 }
01021 
01022 //----------------------------------------------------------------------------
01023 void fei::VectorSpace::getFields(std::vector<int>& fieldIDs)
01024 {
01025   unsigned numFields = fieldDatabase_.size();
01026 
01027   fieldIDs.resize(numFields);
01028 
01029   fei::copyKeysToArray<std::map<int,unsigned> >(fieldDatabase_, numFields,
01030                                                     &fieldIDs[0]);
01031 }
01032 
01033 //----------------------------------------------------------------------------
01034 size_t fei::VectorSpace::getNumIDTypes()
01035 {
01036   return(idTypes_.size());
01037 }
01038 
01039 //----------------------------------------------------------------------------
01040 void fei::VectorSpace::getIDTypes(std::vector<int>& idTypes) const
01041 {
01042   size_t numIDTypes = idTypes_.size();
01043   idTypes.resize(numIDTypes);
01044   for(size_t i=0; i<numIDTypes; ++i) idTypes[i] = idTypes_[i];
01045 }
01046 
01047 //----------------------------------------------------------------------------
01048 int fei::VectorSpace::getNumFields(int idType,
01049                                    int ID)
01050 {
01051   int idindex = fei::binarySearch(idType, idTypes_);
01052   if (idindex < 0) return(0);
01053 
01054   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
01055   if (record == NULL) {
01056     ERReturn(-1);
01057   }
01058 
01059   return( record->getFieldMask()->getNumFields() );
01060 }
01061 
01062 //----------------------------------------------------------------------------
01063 bool fei::VectorSpace::isLocal(int idType, int ID)
01064 {
01065   int idindex = fei::binarySearch(idType, idTypes_);
01066   if (idindex < 0) return(false);
01067 
01068   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
01069   if (record == NULL) {
01070     return(false);
01071   }
01072 
01073   return(true);
01074 }
01075 
01076 //----------------------------------------------------------------------------
01077 bool fei::VectorSpace::isLocallyOwned(int idType, int ID)
01078 {
01079   int idindex = fei::binarySearch(idType, idTypes_);
01080   if (idindex < 0) return(false);
01081 
01082   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
01083   if (record == NULL) {
01084     return(false);
01085   }
01086   if (record->getOwnerProc() == fei::localProc(comm_)) {
01087     return(true);
01088   }
01089 
01090   return(false);
01091 }
01092 
01093 //----------------------------------------------------------------------------
01094 unsigned fei::VectorSpace::getFieldSize(int fieldID)
01095 {
01096   std::map<int,unsigned>::const_iterator
01097     f_iter = fieldDatabase_.find(fieldID);
01098 
01099   if (f_iter == fieldDatabase_.end()) {
01100     FEI_OSTRINGSTREAM osstr;
01101     osstr << "fei::VectorSpace";
01102     if (name_.length() > 0) {
01103       osstr << "(name: "<<name_<<")";
01104     }
01105     osstr << "::getFieldSize: fieldID " << fieldID << " not found.";
01106     throw std::runtime_error(osstr.str());
01107   }
01108 
01109   return((*f_iter).second);
01110 }
01111 
01112 //----------------------------------------------------------------------------
01113 void fei::VectorSpace::getFields(int idType, int ID, std::vector<int>& fieldIDs)
01114 {
01115   int idindex = fei::binarySearch(idType, idTypes_);
01116   if (idindex < 0) {
01117     fieldIDs.resize(0);
01118     return;
01119   }
01120 
01121   fei::Record<int>* record = recordCollections_[idindex]->getRecordWithID(ID);
01122   if (record == NULL) {
01123     voidERReturn;
01124   }
01125 
01126   fei::FieldMask* fieldMask = record->getFieldMask();
01127   std::vector<int>& maskFieldIDs = fieldMask->getFieldIDs();
01128   int numFields = maskFieldIDs.size();
01129   fieldIDs.resize(numFields);
01130   for(int i=0; i<numFields; ++i) {
01131     fieldIDs[i] = maskFieldIDs[i];
01132   }
01133 }
01134 
01135 //----------------------------------------------------------------------------
01136 void fei::VectorSpace::getGlobalIndexOffsets(std::vector<int>& globalOffsets) const
01137 {
01138   globalOffsets = globalOffsets_;
01139 }
01140 
01141 //----------------------------------------------------------------------------
01142 void fei::VectorSpace::getGlobalBlkIndexOffsets(std::vector<int>& globalBlkOffsets) const
01143 {
01144   globalBlkOffsets = globalIDOffsets_;
01145 }
01146 
01147 //----------------------------------------------------------------------------
01148 int fei::VectorSpace::getOwnerProcPtIndex(int globalIndex)
01149 {
01150   if (globalIndex < 0) return(-1);
01151 
01152   unsigned len = globalOffsets_.size();
01153   for(int i=0; i<(int)(len-1); ++i) {
01154     if (globalIndex < globalOffsets_[i+1]) {
01155       return(i);
01156     }
01157   }
01158 
01159   return(-1);
01160 }
01161 
01162 //----------------------------------------------------------------------------
01163 int fei::VectorSpace::getOwnerProcBlkIndex(int globalIndex)
01164 {
01165   if (globalIndex < 0) return(-1);
01166 
01167   unsigned len = globalOffsets_.size();
01168   for(int i=0; i<(int)(len-1); ++i) {
01169     if (globalIndex < globalIDOffsets_[i+1]) {
01170       return(i);
01171     }
01172   }
01173 
01174   return(-1);
01175 }
01176 
01177 //----------------------------------------------------------------------------
01178 int fei::VectorSpace::getNumOwnedAndSharedIDs(int idType)
01179 {
01180   int idx = fei::binarySearch(idType, idTypes_);
01181   if (idx < 0) return(0);
01182 
01183   return( recordCollections_[idx]->getNumRecords() );
01184 }
01185 
01186 //----------------------------------------------------------------------------
01187 int fei::VectorSpace::getNumOwnedIDs(int idType)
01188 {
01189   int idx = fei::binarySearch(idType, idTypes_);
01190   if (idx < 0) return(0);
01191 
01192   fei::RecordAttributeCounter attrCounter(fei::localProc(comm_));
01193   runRecords(attrCounter);
01194 
01195   return( attrCounter.numLocallyOwnedIDs_ );
01196 }
01197 
01198 //----------------------------------------------------------------------------
01199 int fei::VectorSpace::getNumSharedIDs(int idType, int& numShared)
01200 {
01201   numShared = sharedIDTables_[idType].getSharedIDs().size();
01202 
01203   return(0);
01204 }
01205 
01206 //----------------------------------------------------------------------------
01207 int fei::VectorSpace::getOwnedAndSharedIDs(int idType,
01208                                       int lenList,
01209                                       int* IDs,
01210                                       int& numLocalIDs)
01211 {
01212   int idx = fei::binarySearch(idType, idTypes_);
01213   if (idx < 0) return(-1);
01214 
01215   snl_fei::RecordCollection* records = recordCollections_[idx];
01216 
01217   snl_fei::RecordCollection::map_type& rmap = records->getRecords();
01218 
01219   numLocalIDs = rmap.size();
01220 
01221   fei::copyKeysToArray(rmap, lenList, IDs);
01222 
01223   return(0);
01224 }
01225 
01226 //----------------------------------------------------------------------------
01227 int fei::VectorSpace::getOwnedIDs(int idType,
01228                                          int lenList,
01229                                          int* IDs,
01230                                          int& numLocalIDs)
01231 {
01232   int idx = fei::binarySearch(idType, idTypes_);
01233   if (idx < 0) return(-1);
01234 
01235   snl_fei::RecordCollection* records = recordCollections_[idx];
01236 
01237   snl_fei::RecordCollection::map_type& rmap = records->getRecords();
01238 
01239   numLocalIDs = 0;
01240 
01241   snl_fei::RecordCollection::map_type::iterator
01242     r_iter = rmap.begin(),
01243     r_end  = rmap.end();
01244 
01245   for(; r_iter != r_end; ++r_iter) {
01246     fei::Record<int>& thisrecord = *((*r_iter).second);
01247 
01248     if (thisrecord.getOwnerProc() == fei::localProc(comm_)) {
01249       if (numLocalIDs < lenList) {
01250         IDs[numLocalIDs] = thisrecord.getID();
01251       }
01252       ++numLocalIDs;
01253     }
01254   }
01255 
01256   return(0);
01257 }
01258 
01259 //----------------------------------------------------------------------------
01260 int fei::VectorSpace::getNumIndices_SharedAndOwned() const
01261 {
01262   return(eqnNumbers_.size());
01263 }
01264 
01265 //----------------------------------------------------------------------------
01266 int fei::VectorSpace::getIndices_SharedAndOwned(std::vector<int>& globalIndices) const
01267 {
01268   if (eqnNumbers_.size() == 0) {
01269     globalIndices.resize(0);
01270     return(0);
01271   }
01272 
01273   size_t numIndices = eqnNumbers_.size();
01274   const int* indicesPtr = &eqnNumbers_[0];
01275 
01276   globalIndices.resize(numIndices);
01277   for(size_t i=0; i<numIndices; ++i) {
01278     globalIndices[i] = indicesPtr[i];
01279   }
01280 
01281   return(0);
01282 }
01283 
01284 //----------------------------------------------------------------------------
01285 int fei::VectorSpace::getNumBlkIndices_SharedAndOwned(int& numBlkIndices) const
01286 {
01287   numBlkIndices = 0;
01288   for(size_t i=0; i<recordCollections_.size(); ++i) {
01289     numBlkIndices += recordCollections_[i]->getNumRecords();
01290   }
01291 
01292   return(0);
01293 }
01294 
01295 //----------------------------------------------------------------------------
01296 int fei::VectorSpace::getBlkIndices_SharedAndOwned(int lenBlkIndices,
01297                                                    int* globalBlkIndices,
01298                                                    int* blkSizes,
01299                                                    int& numBlkIndices)
01300 {
01301   if (!sharedRecordsSynchronized_) {
01302     numBlkIndices = 0;
01303     return(-1);
01304   }
01305 
01306   fei::BlkIndexAccessor blkIndAccessor(lenBlkIndices,
01307                                   globalBlkIndices, blkSizes);
01308   runRecords(blkIndAccessor);
01309 
01310   numBlkIndices = blkIndAccessor.numBlkIndices_;
01311 
01312   return(0);
01313 }
01314 
01315 //----------------------------------------------------------------------------
01316 int fei::VectorSpace::getGlobalNumIndices() const
01317 {
01318   if (globalOffsets_.size() < 1) return(0);
01319   return(globalOffsets_[globalOffsets_.size()-1]);
01320 }
01321 
01322 //----------------------------------------------------------------------------
01323 int fei::VectorSpace::getNumIndices_Owned() const
01324 {
01325   if (!sharedRecordsSynchronized_) {
01326     return(-1);
01327   }
01328 
01329   int localProc = fei::localProc(comm_);
01330   int numIndices = globalOffsets_[localProc+1]-globalOffsets_[localProc];
01331 
01332   return(numIndices);
01333 }
01334 
01335 //----------------------------------------------------------------------------
01336 int fei::VectorSpace::getIndices_Owned(std::vector<int>& globalIndices) const
01337 {
01338   if (!sharedRecordsSynchronized_) {
01339     globalIndices.clear();
01340     return(-1);
01341   }
01342 
01343   int localProc = fei::localProc(comm_);
01344   int numIndices = globalOffsets_[localProc+1]-globalOffsets_[localProc];
01345 
01346   globalIndices.resize(numIndices);
01347 
01348   int firstOffset = globalOffsets_[localProc];
01349   for(int i=0; i<numIndices; ++i) {
01350     globalIndices[i] = firstOffset+i;
01351   }
01352 
01353   return(0);
01354 }
01355 
01356 //----------------------------------------------------------------------------
01357 int fei::VectorSpace::getIndices_Owned(int lenIndices, int* globalIndices, int& numIndices) const
01358 {
01359   if (!sharedRecordsSynchronized_) {
01360     numIndices = 0;
01361     return(-1);
01362   }
01363 
01364   int localProc = fei::localProc(comm_);
01365   numIndices = globalOffsets_[localProc+1]-globalOffsets_[localProc];
01366 
01367   int len = lenIndices >= numIndices ? numIndices : lenIndices;
01368 
01369   int firstOffset = globalOffsets_[localProc];
01370   for(int i=0; i<len; ++i) {
01371     globalIndices[i] = firstOffset+i;
01372   }
01373 
01374   return(0);
01375 }
01376 
01377 //----------------------------------------------------------------------------
01378 int fei::VectorSpace::getNumBlkIndices_Owned() const
01379 {
01380   if (!sharedRecordsSynchronized_) {
01381     return(-1);
01382   }
01383 
01384   int localProc = fei::localProc(comm_);
01385   int numBlkIndices = globalIDOffsets_[localProc+1]-globalIDOffsets_[localProc];
01386 
01387   return(numBlkIndices);
01388 }
01389 
01390 //----------------------------------------------------------------------------
01391 int fei::VectorSpace::getGlobalNumBlkIndices() const
01392 {
01393   int numBlkIndices = 0;
01394   unsigned len = globalIDOffsets_.size();
01395   if (len > 0) {
01396     numBlkIndices = globalIDOffsets_[len-1];
01397   }
01398 
01399   return(numBlkIndices);
01400 }
01401 
01402 //----------------------------------------------------------------------------
01403 int fei::VectorSpace::getBlkIndices_Owned(int lenBlkIndices,
01404                                           int* globalBlkIndices,
01405                                           int* blkSizes,
01406                                           int& numBlkIndices)
01407 {
01408   if (!sharedRecordsSynchronized_) {
01409     numBlkIndices = 0;
01410     return(-1);
01411   }
01412 
01413   int localProc = fei::localProc(comm_);
01414   fei::BlkIndexAccessor blkIndAccessor(localProc, lenBlkIndices,
01415                                            globalBlkIndices, blkSizes);
01416   runRecords(blkIndAccessor);
01417 
01418   numBlkIndices = blkIndAccessor.numBlkIndices_;
01419 
01420   return(0);
01421 }
01422 
01423 //----------------------------------------------------------------------------
01424 int fei::VectorSpace::getRecordCollection(int idType,
01425                                           snl_fei::RecordCollection*& records)
01426 {
01427   int idx = fei::binarySearch(idType, idTypes_);
01428   if (idx < 0) return(-1);
01429 
01430   records = recordCollections_[idx];
01431   return(0);
01432 }
01433 
01434 //----------------------------------------------------------------------------
01435 int fei::VectorSpace::getRecordCollection(int idType,
01436                                           const snl_fei::RecordCollection*& records) const
01437 {
01438   int idx = fei::binarySearch(idType, idTypes_);
01439   if (idx < 0) return(-1);
01440 
01441   records = recordCollections_[idx];
01442   return(0);
01443 }
01444 
01445 //----------------------------------------------------------------------------
01446 void fei::VectorSpace::setOwners_lowestSharing()
01447 {
01448   //first, add localProc to each of the sharing-proc lists, in case it wasn't
01449   //included via initSharedIDs().
01450   int localProc = fei::localProc(comm_);
01451 
01452   std::map<int, fei::SharedIDs<int> >::iterator
01453     s_iter = sharedIDTables_.begin(), s_end = sharedIDTables_.end();
01454 
01455   for(; s_iter != s_end; ++s_iter) {
01456     fei::SharedIDs<int>::map_type& shid_table = s_iter->second.getSharedIDs();
01457 
01458     fei::SharedIDs<int>::map_type::iterator
01459       t_iter = shid_table.begin(),
01460       t_end = shid_table.end();
01461 
01462     for(; t_iter != t_end; ++t_iter) {
01463       std::set<int>& shProcs = t_iter->second;
01464       shProcs.insert(localProc);
01465     }
01466   }
01467 
01468   //now set the owningProcs for the SharedIDs records, and the owning procs on
01469   //the appropriate records in the recordCollections. Set the owner to be the
01470   //lowest-numbered sharing proc in all cases.
01471   s_iter = sharedIDTables_.begin();
01472   for(; s_iter != s_end; ++s_iter) {
01473     fei::SharedIDs<int>::map_type& shid_table = s_iter->second.getSharedIDs();
01474 
01475     std::vector<int>& owningProcs = s_iter->second.getOwningProcs();
01476 
01477     int len = shid_table.size();
01478     owningProcs.resize(len);
01479     int j = 0;
01480 
01481     fei::SharedIDs<int>::map_type::iterator
01482       t_iter = shid_table.begin(),
01483       t_end = shid_table.end();
01484 
01485     for(; t_iter != t_end; ++t_iter, ++j) {
01486       std::set<int>& shProcs = (*t_iter).second;
01487       int lowest = *(shProcs.begin());
01488       owningProcs[j] = lowest;
01489     }
01490     
01491     int idx = fei::binarySearch(s_iter->first, idTypes_);
01492     if (idx < 0) {
01493       throw std::runtime_error("internal error in fei::VectorSapce::setOwners_lowestSharing");
01494     }
01495 
01496     if (output_level_ >= fei::FULL_LOGS && output_stream_ != NULL) {
01497       recordCollections_[idx]->setDebugOutput(output_stream_);
01498     }
01499 
01500     recordCollections_[idx]->setOwners_lowestSharing(s_iter->second);
01501   }
01502 }
01503 
01504 //----------------------------------------------------------------------------
01505 int fei::VectorSpace::calculateGlobalIndices()
01506 {
01507   int localProc = fei::localProc(comm_);
01508   int numProcs = fei::numProcs(comm_);
01509   std::vector<int> localOffsets(numProcs +1, 0);
01510   std::vector<int> localIDOffsets(numProcs+1, 0);
01511   globalOffsets_.resize(numProcs + 1, 0);
01512 
01513   globalIDOffsets_.assign(globalIDOffsets_.size(), 0);
01514 
01515   //first we'll calculate the number of local degrees of freedom, and the
01516   //number of local identifiers.
01517   fei::RecordAttributeCounter counter(localProc);
01518   runRecords(counter);
01519 
01520   int numLocalDOF = counter.numLocalDOF_;
01521   int numLocalIDs = counter.numLocalIDs_;
01522   int numRemoteSharedDOF = counter.numRemoteSharedDOF_;
01523 
01524   eqnNumbers_.resize(numLocalDOF+numRemoteSharedDOF);
01525 
01526   localOffsets[localProc] = numLocalDOF;
01527   CHK_MPI( fei::GlobalMax(comm_, localOffsets, globalOffsets_) );
01528 
01529   localIDOffsets[localProc] = numLocalIDs;
01530   CHK_MPI( fei::GlobalMax(comm_, localIDOffsets, globalIDOffsets_) );
01531 
01532   //Now globalOffsets_ contains numLocalDOF for proc i, in the i-th position.
01533   //(and similarly for globalIDOffsets_)
01534   //So all we need to do is turn that data into global-offsets (i.e., the
01535   //starting global-offset for each processor).
01536   int localOffset = 0;
01537   int localIDOffset = 0;
01538   for(int p=0; p<numProcs; ++p) {
01539     numLocalDOF = globalOffsets_[p];
01540     globalOffsets_[p] = localOffset;
01541     localOffset += numLocalDOF;
01542     numLocalIDs = globalIDOffsets_[p];
01543     globalIDOffsets_[p] = localIDOffset;
01544     localIDOffset += numLocalIDs;
01545   }
01546   globalOffsets_[numProcs] = localOffset;
01547   globalIDOffsets_[numProcs] = localIDOffset;
01548 
01549   firstLocalOffset_ = globalOffsets_[localProc];
01550   lastLocalOffset_ = globalOffsets_[localProc+1] - 1;
01551 
01552   if (output_level_ >= fei::BRIEF_LOGS && output_stream_ != NULL) {
01553     FEI_OSTREAM& os = *output_stream_;
01554     os <<dbgprefix_<< "  firstLocalOffset_: " << firstLocalOffset_ << ", lastLocalOffset_: "
01555       << lastLocalOffset_ << FEI_ENDL;
01556   }
01557 
01558   //Now we're ready to set the equation-numbers on all local records.
01559   CHK_ERR( setLocalEqnNumbers() );
01560 
01561   return(0);
01562 }
01563 
01564 //----------------------------------------------------------------------------
01565 int fei::VectorSpace::synchronizeSharedRecords()
01566 {
01567   if (output_level_ >= fei::FULL_LOGS && output_stream_ != NULL) {
01568     FEI_OSTREAM& os = *output_stream_;
01569     os<<dbgprefix_<<"#synchronizeSharedRecords num-field-masks: "<<fieldMasks_.size()<<FEI_ENDL;
01570     for(unsigned fm=0; fm<fieldMasks_.size(); ++fm) {
01571       os << dbgprefix_<<"#     maskID["<<fm<<"]: " << fieldMasks_[fm]->getMaskID() << FEI_ENDL;
01572     }
01573   }
01574 
01575   bool safetyCheck = checkSharedIDs_;
01576 
01577   int numProcs = fei::numProcs(comm_);
01578   int localProc = fei::localProc(comm_);
01579   int numShTables = sharedIDTables_.size();
01580 
01581   if (numProcs < 2) {
01582     sharedRecordsSynchronized_ = true;
01583     return(0);
01584   }
01585 
01586   if (safetyCheck) {
01587     if (localProc == 0 && output_level_ >= fei::BRIEF_LOGS) {
01588       FEI_COUT << "fei::VectorSpace: global consistency-check of shared ID"
01589            << " data (involves all-to-all communication). This is done "
01590            << "only if requested by parameter 'FEI_CHECK_SHARED_IDS true'."
01591            << FEI_ENDL;
01592     }
01593 
01594     int totalNumShTables = 0;
01595     CHK_ERR( fei::GlobalSum(comm_, numShTables, totalNumShTables) );
01596     if (totalNumShTables != numShTables * numProcs) {
01597       //not all processors have the same number of shared-id tables. This means
01598       //that one or more processors is 'empty', which may not be an error. But
01599       //it does mean that the safety check can't be performed because the safety
01600       //check involves all-to-all communication and if one or more processors
01601       //don't participate then we'll hang...
01602       safetyCheck = false;
01603     }
01604   }
01605 
01606   //Create a list of fei::comm_maps which will be the communication-patterns
01607   //for each of the sharedIDTables. A sharedIDTable maps lists of processors
01608   //to each shared ID. The communication-pattern will be a transpose of
01609   //that, mapping lists of IDs to owning or sharing processors.
01610   //
01611 
01612   int local_err = 0;
01613 
01614   for(size_t i=0; i<idTypes_.size(); ++i) {
01615 
01616     fei::SharedIDs<int>& shared = getSharedIDs_private(idTypes_[i]);
01617     
01618     //now create/initialize ownerPatterns which map owning processors to lists
01619     //of ids that are shared locally, and sharerPatterns which map sharing
01620     //processors to lists of ids that are owned locally.
01621     fei::comm_map* ownerPattern = new fei::comm_map(1, numProcs);
01622 
01623     ownerPatterns_[idTypes_[i]] = ownerPattern;
01624 
01625     fei::comm_map* sharerPattern = 
01626             new fei::comm_map(1, numProcs);
01627 
01628     sharerPatterns_[idTypes_[i]] = sharerPattern;
01629 
01630     std::vector<int>& owningProcs = shared.getOwningProcs();
01631 
01632     fei::SharedIDs<int>::map_type& shtable = shared.getSharedIDs();
01633 
01634     fei::SharedIDs<int>::map_type::iterator
01635       s_iter = shtable.begin(),
01636       s_end = shtable.end();
01637     
01638     int j = 0;
01639     for(; s_iter != s_end; ++s_iter, ++j) {
01640       int ID = (*s_iter).first;
01641       std::set<int>& shProcs = (*s_iter).second;
01642 
01643       int owner = owningProcs[j];
01644 
01645       if (owner == localProc) {
01646         std::set<int>::const_iterator
01647           p_iter = shProcs.begin(),
01648           p_end = shProcs.end();
01649         for(; p_iter != p_end; ++p_iter) {
01650           if (*p_iter != localProc) {
01651             sharerPattern->addIndices(*p_iter, 1, &ID);
01652           }
01653         }
01654       }
01655       else {
01656         ownerPattern->addIndices(owner, 1, &ID);
01657       }
01658     }
01659 
01660     if (safetyCheck) {
01661       fei::comm_map* checkPattern = NULL;
01662       CHK_ERR( fei::mirrorCommPattern(comm_, ownerPattern, checkPattern) );
01663       fei::Barrier(comm_);
01664       if (output_level_ >= fei::FULL_LOGS && output_stream_ != NULL) {
01665         FEI_OSTREAM& os = *output_stream_;
01666         fei::comm_map::map_type& owner_map = ownerPattern->getMap();
01667         int numKeys = owner_map.size();
01668         fei::comm_map::map_type::const_iterator
01669           omap_iter = owner_map.begin(),
01670           omap_end = owner_map.end();
01671 
01672         os << dbgprefix_<<"#  synchronizeSharedRecords" << FEI_ENDL
01673            << dbgprefix_<<"#  ownerPattern, num-procs-to-send-to: " << numKeys << FEI_ENDL;
01674         for(int sk=0; omap_iter != omap_end; ++sk, ++omap_iter) {
01675           os << dbgprefix_<<"#    sendProc["<<sk<<"]: " << omap_iter->first << " IDs: ";
01676           fei::comm_map::row_type::const_iterator
01677             val_iter = omap_iter->second->begin(),
01678             val_end =  omap_iter->second->end();
01679           for(; val_iter != val_end; ++val_iter) {
01680             os << *val_iter << " ";
01681           }
01682           os << FEI_ENDL;
01683         }
01684 
01685         fei::comm_map::map_type& check_map = checkPattern->getMap();
01686         int numCheckKeys = check_map.size();
01687         fei::comm_map::map_type::const_iterator
01688           cmap_iter = check_map.begin(),
01689           cmap_end = check_map.end();
01690 
01691         os <<dbgprefix_<< "#  synchronizeSharedRecords" << FEI_ENDL
01692            <<dbgprefix_<< "#  checkPattern (send mirror), num-procs: "
01693            << numCheckKeys << FEI_ENDL;
01694         for(int sk=0; cmap_iter != cmap_end; ++sk, ++cmap_iter) {
01695           os <<dbgprefix_<< "#    proc["<<sk<<"]: " << cmap_iter->first << " IDs: ";
01696           fei::comm_map::row_type::const_iterator
01697           val_iter = cmap_iter->second->begin(),
01698             val_end = cmap_iter->second->end();
01699           for(; val_iter != val_end; ++val_iter) {
01700             os << *val_iter << " ";
01701           }
01702           os << FEI_ENDL;
01703         }
01704       }
01705 
01706       int err = 0;
01707       bool quiet = false;
01708       if (!checkPattern->equal(*sharerPattern, quiet)) {
01709         //FEI_CERR << "shared-ID safety-check failed..."
01710         //     << FEI_ENDL;
01711         err = -1;
01712       }
01713       int globalErr = 0;
01714       CHK_ERR( fei::GlobalSum(comm_, err, globalErr) );
01715 
01716       delete checkPattern;
01717 
01718       if (globalErr != 0) {
01719         return(globalErr);
01720       }
01721     }
01722 
01723     local_err += exchangeFieldInfo(ownerPattern, sharerPattern,
01724                                    recordCollections_[i], fieldMasks_);
01725   }
01726 
01727   int global_err = 0;
01728   CHK_ERR( fei::GlobalSum(comm_, local_err, global_err) );
01729 
01730   if (global_err != 0) {
01731     ERReturn(-1);
01732   }
01733 
01734   sharedRecordsSynchronized_ = true;
01735 
01736   return(0);
01737 }
01738 
01739 //----------------------------------------------------------------------------
01740 int fei::VectorSpace::exchangeGlobalIndices()
01741 {
01742   for(size_t i=0; i<idTypes_.size(); ++i) {
01743     snl_fei::RecordMsgHandler recmsghndlr(fei::localProc(comm_),
01744                                  recordCollections_[i], *ptBlkMap_,
01745                                           fieldMasks_, eqnNumbers_);
01746     recmsghndlr.setTask(snl_fei::RecordMsgHandler::_EqnNumbers_);
01747 
01748     recmsghndlr.setSendPattern(sharerPatterns_[idTypes_[i]]);
01749     recmsghndlr.setRecvPattern(ownerPatterns_[idTypes_[i]]);
01750     CHK_ERR( fei::exchange(comm_, &recmsghndlr) );
01751   }
01752 
01753   return(0);
01754 }
01755 
01756 //----------------------------------------------------------------------------
01757 void fei::VectorSpace::runRecords(fei::Record_Operator<int>& record_op)
01758 {
01759   for(size_t rec=0; rec<recordCollections_.size(); ++rec) {
01760     snl_fei::RecordCollection* records = recordCollections_[rec];
01761 
01762     snl_fei::RecordCollection::map_type& rmap = records->getRecords();
01763 
01764     snl_fei::RecordCollection::map_type::iterator
01765       r_iter = rmap.begin(),
01766       r_end  = rmap.end();
01767 
01768     for(; r_iter != r_end; ++r_iter) {
01769       fei::Record<int>& thisrecord = *((*r_iter).second);
01770 
01771       record_op(thisrecord);
01772     }
01773   }
01774 }
01775 
01776 //----------------------------------------------------------------------------
01777 int fei::VectorSpace::setLocalEqnNumbers()
01778 {
01779   int proc = fei::localProc(comm_);
01780   int eqnNumber = firstLocalOffset_;
01781   int idNumber = globalIDOffsets_[proc];
01782 
01783   int numProcs = fei::numProcs(comm_);
01784   int localProc = fei::localProc(comm_);
01785 
01786   if (ptBlkMap_ != NULL) {
01787     delete ptBlkMap_;
01788   }
01789   ptBlkMap_ = new snl_fei::PointBlockMap;
01790 
01791   int maxNumIndices = 0;
01792   for(unsigned i=0; i<fieldMasks_.size(); ++i) {
01793     if (fieldMasks_[i]->getNumIndices() > maxNumIndices) {
01794       maxNumIndices = fieldMasks_[i]->getNumIndices();
01795     }
01796   }
01797 
01798   if (maxNumIndices == 1) {
01799     ptBlkMap_->setPtEqualBlk();
01800   }
01801 
01802   FEI_OSTREAM* id2eqnStream = NULL;
01803   if (output_level_ >= fei::BRIEF_LOGS) {
01804     std::string path = fei::LogManager::getLogManager().getOutputPath();
01805     if (path == "") path = ".";
01806     FEI_OSTRINGSTREAM osstr;
01807     osstr << path << "/fei_ID2Eqn";
01808     if (name_.size() > 0) osstr << "_" << name_;
01809     osstr << "." <<numProcs<< "." << localProc;
01810 
01811     id2eqnStream = new FEI_OFSTREAM(osstr.str().c_str(), IOS_OUT);
01812     FEI_OSTREAM& os = *id2eqnStream;
01813     os << "# Each line contains:\n#   ID   blk-eqn   eqn" << FEI_ENDL;
01814   }
01815 
01816   int eqnNumberOffset = 0;
01817   int maxNumDOF = 0;
01818 
01819   for(size_t rec=0; rec<recordCollections_.size(); ++rec) {
01820     snl_fei::RecordCollection* records = recordCollections_[rec];
01821 
01822     snl_fei::RecordCollection::map_type& rmap = records->getRecords();
01823 
01824     snl_fei::RecordCollection::map_type::iterator
01825       r_iter = rmap.begin(),
01826       r_end = rmap.end();
01827 
01828     int* eqnNumPtr = &eqnNumbers_[0];
01829 
01830     for(; r_iter != r_end; ++r_iter) {
01831       fei::Record<int>& thisrecord = *((*r_iter).second);
01832 
01833       fei::FieldMask* mask = thisrecord.getFieldMask();
01834 
01835       thisrecord.setOffsetIntoEqnNumbers(eqnNumberOffset);
01836 
01837       int owner = thisrecord.getOwnerProc();
01838       if (owner == proc) {
01839         thisrecord.setNumber(idNumber++);
01840       }
01841 
01842       int* eqnNumbers = eqnNumPtr
01843                       + thisrecord.getOffsetIntoEqnNumbers();
01844 
01845       int numDOF = mask->getNumIndices();
01846       eqnNumberOffset += numDOF;
01847 
01848       if (output_level_ >= fei::BRIEF_LOGS) {
01849         for(int ii=0; ii<numDOF; ++ii) {
01850           if (isLogEqn(eqnNumber+ii) && output_stream_ != NULL) {
01851             FEI_OSTREAM& os = *output_stream_;
01852             os <<dbgprefix_<< "setLocalEqnNumbers: ID "<<thisrecord.getID()
01853                << " <--> eqn " << eqnNumber+ii<<FEI_ENDL;
01854           }
01855         }
01856       }
01857 
01858       if (owner == proc) {
01859         int offset = 0;
01860         for(int n=0; n<numDOF; ++n) {
01861           eqnNumbers[offset++] = eqnNumber++;
01862         }
01863       }
01864 
01865       if (numDOF > maxNumDOF) maxNumDOF = numDOF;
01866 
01867       if (owner == proc) {
01868         int thiseqn = eqnNumber-numDOF;
01869         int thisrecordnumber = thisrecord.getNumber();
01870         if (maxNumIndices > 1) {
01871           CHK_ERR( ptBlkMap_->setEqn(thiseqn, thisrecordnumber, numDOF) );
01872           if (numDOF > 1) {
01873             for(int i=1; i<numDOF; ++i) {
01874               CHK_ERR( ptBlkMap_->setEqn(thiseqn+i, thisrecordnumber, numDOF) );
01875             }
01876           }
01877         }
01878       }
01879 
01880       if (id2eqnStream != NULL) {
01881         if (owner == proc) {
01882           FEI_OSTREAM& os = *id2eqnStream;
01883           for(int n=0; n<numDOF; ++n) {
01884             os << thisrecord.getID() << " " << thisrecord.getNumber() << " "
01885                << eqnNumber-numDOF+n<<FEI_ENDL;
01886           }
01887         }
01888       }
01889     }
01890 
01891   }
01892 
01893   ptBlkMap_->setMaxBlkEqnSize(maxNumDOF);
01894 
01895   int globalMaxNumDOF;
01896   CHK_ERR( fei::GlobalMax(comm_, maxNumDOF, globalMaxNumDOF) );
01897 
01898   if (globalMaxNumDOF == 1) {
01899     ptBlkMap_->setPtEqualBlk();
01900   }
01901 
01902   delete id2eqnStream;
01903 
01904   return(0);
01905 }
01906 
01907 //----------------------------------------------------------------------------
01908 int fei::VectorSpace::exchangeFieldInfo(fei::comm_map* ownerPattern,
01909                                         fei::comm_map* sharerPattern,
01910                                   snl_fei::RecordCollection* recordCollection,
01911                                       std::vector<fei::FieldMask*>& fieldMasks)
01912 {
01913   //ownerPattern maps owning processors to lists of IDs that we (the local
01914   //processor) share.
01915   //
01916   //sharerPattern maps sharing processors to lists of IDs that we own.
01917   //
01918   //In this function we need to perform these tasks:
01919   //1. processors exchange and combine their sets of field-masks so that all
01920   //   processors will have the super-set of field-masks that they need.
01921   //2. sharing processors send maskIDs for their shared IDs to the owning
01922   //   processors. The owning processors then combine masks as necessary to
01923   //   make sure each shared ID has the union of field-masks that are held by
01924   //   each of the sharing processors. all owning processors now send their
01925   //   maskIDs out to the sharing processors. At this point all processors
01926   //   should have the same view of the masks that belong on shared IDs.
01927   //3. exchange info describing inactive DOF offsets for shared records.
01928   //
01929 
01930   int numProcs = fei::numProcs(comm_);
01931   if (numProcs < 2) return(0);
01932 
01933   if (ptBlkMap_ == 0) ptBlkMap_ = new snl_fei::PointBlockMap;
01934 
01935   snl_fei::RecordMsgHandler recMsgHandler(fei::localProc(comm_),
01936                                           recordCollection,
01937                                           *ptBlkMap_,
01938                                           fieldMasks, eqnNumbers_);
01939 
01940   //Step 1a.
01941   recMsgHandler.setTask(snl_fei::RecordMsgHandler::_FieldMasks_);
01942 
01943   recMsgHandler.setSendPattern(ownerPattern);
01944   recMsgHandler.setRecvPattern(sharerPattern);
01945   CHK_ERR( fei::exchange(comm_, &recMsgHandler) );
01946 
01947   //Step 2a.
01948   recMsgHandler.setTask(snl_fei::RecordMsgHandler::_MaskIDs_);
01949 
01950   recMsgHandler.setSendPattern(ownerPattern);
01951   recMsgHandler.setRecvPattern(sharerPattern);
01952   CHK_ERR( fei::exchange(comm_, &recMsgHandler) );
01953 
01954   //Step 1b.
01955   recMsgHandler.setTask(snl_fei::RecordMsgHandler::_FieldMasks_);
01956 
01957   recMsgHandler.setSendPattern(sharerPattern);
01958   recMsgHandler.setRecvPattern(ownerPattern);
01959   CHK_ERR( fei::exchange(comm_, &recMsgHandler) );
01960 
01961   //Step 2b.
01962   recMsgHandler.setTask(snl_fei::RecordMsgHandler::_MaskIDs_);
01963 
01964   recMsgHandler.setSendPattern(sharerPattern);
01965   recMsgHandler.setRecvPattern(ownerPattern);
01966   CHK_ERR( fei::exchange(comm_, &recMsgHandler) );
01967 
01968   return(0);
01969 }
01970 
01971 //----------------------------------------------------------------------------
01972 fei::FieldDofMap<int>&
01973 fei::VectorSpace::getFieldDofMap()
01974 {
01975   return fieldDofMap_;
01976 }
01977 
01978 //----------------------------------------------------------------------------
01979 void fei::VectorSpace::setName(const char* name)
01980 {
01981   if (name == NULL) return;
01982 
01983   if (name_ == name) return;
01984 
01985   name_ = name;
01986   dbgprefix_ = "VecSpc_"+name_+": ";
01987 }
01988 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends