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

Generated on Tue Jul 13 09:27:46 2010 for FEI by  doxygen 1.4.7