Epetra Package Browser (Single Doxygen Collection) Development
Epetra_FEVector.cpp
Go to the documentation of this file.
00001 
00002 //@HEADER
00003 // ************************************************************************
00004 //
00005 //               Epetra: Linear Algebra Services Package
00006 //                 Copyright 2011 Sandia Corporation
00007 //
00008 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00009 // the U.S. Government retains certain rights in this software.
00010 //
00011 // Redistribution and use in source and binary forms, with or without
00012 // modification, are permitted provided that the following conditions are
00013 // met:
00014 //
00015 // 1. Redistributions of source code must retain the above copyright
00016 // notice, this list of conditions and the following disclaimer.
00017 //
00018 // 2. Redistributions in binary form must reproduce the above copyright
00019 // notice, this list of conditions and the following disclaimer in the
00020 // documentation and/or other materials provided with the distribution.
00021 //
00022 // 3. Neither the name of the Corporation nor the names of the
00023 // contributors may be used to endorse or promote products derived from
00024 // this software without specific prior written permission.
00025 //
00026 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00027 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00028 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00029 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00030 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00031 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00032 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00033 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00034 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00035 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00036 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00037 //
00038 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00039 //
00040 // ************************************************************************
00041 //@HEADER
00042 
00043 #include <Epetra_ConfigDefs.h>
00044 #include <Epetra_FEVector.h>
00045 
00046 #include <Epetra_LocalMap.h>
00047 #include <Epetra_Comm.h>
00048 #include <Epetra_Map.h>
00049 #include <Epetra_Import.h>
00050 #include <Epetra_Export.h>
00051 #include <Epetra_Util.h>
00052 #include <Epetra_IntSerialDenseVector.h>
00053 #include <Epetra_SerialDenseVector.h>
00054 
00055 #include <algorithm>
00056 
00057 //----------------------------------------------------------------------------
00058 Epetra_FEVector::Epetra_FEVector(const Epetra_BlockMap& map,
00059                                  int numVectors,
00060          bool ignoreNonLocalEntries)
00061   : Epetra_MultiVector(map, numVectors),
00062     myFirstID_(0),
00063     myNumIDs_(0),
00064 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00065     nonlocalIDs_int_(),
00066 #endif
00067 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00068     nonlocalIDs_LL_(),
00069 #endif
00070     nonlocalElementSize_(),
00071     nonlocalCoefs_(),
00072     nonlocalMap_(0),
00073     exporter_(0),
00074     nonlocalVector_(0),
00075     ignoreNonLocalEntries_(ignoreNonLocalEntries)
00076 {
00077   myFirstID_ = map.MinMyGID64();
00078   myNumIDs_ = map.NumMyElements();
00079   nonlocalCoefs_.resize(numVectors);
00080 }
00081 
00082 //----------------------------------------------------------------------------
00083 Epetra_FEVector::Epetra_FEVector(Epetra_DataAccess CV, const Epetra_BlockMap& theMap,
00084                                  double *A, int MyLDA, int theNumVectors,
00085                                  bool ignoreNonLocalEntries)
00086  : Epetra_MultiVector(CV, theMap, A, MyLDA, theNumVectors),
00087     myFirstID_(0),
00088     myNumIDs_(0),
00089 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00090     nonlocalIDs_int_(),
00091 #endif
00092 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00093     nonlocalIDs_LL_(),
00094 #endif
00095     nonlocalElementSize_(),
00096     nonlocalCoefs_(),
00097     nonlocalMap_(0),
00098     exporter_(0),
00099     nonlocalVector_(0),
00100     ignoreNonLocalEntries_(ignoreNonLocalEntries)
00101 {
00102   myFirstID_ = theMap.MinMyGID64();
00103   myNumIDs_ = theMap.NumMyElements();
00104   nonlocalCoefs_.resize(theNumVectors);
00105 }
00106 
00107 //----------------------------------------------------------------------------
00108 Epetra_FEVector::Epetra_FEVector(Epetra_DataAccess CV, const Epetra_BlockMap& theMap,
00109                                  double **ArrayOfPointers, int theNumVectors,
00110                                  bool ignoreNonLocalEntries)
00111  : Epetra_MultiVector(CV, theMap, ArrayOfPointers, theNumVectors),
00112     myFirstID_(0),
00113     myNumIDs_(0),
00114 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00115     nonlocalIDs_int_(),
00116 #endif
00117 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00118     nonlocalIDs_LL_(),
00119 #endif
00120     nonlocalElementSize_(),
00121     nonlocalCoefs_(),
00122     nonlocalMap_(0),
00123     exporter_(0),
00124     nonlocalVector_(0),
00125     ignoreNonLocalEntries_(ignoreNonLocalEntries)
00126 {
00127   myFirstID_ = theMap.MinMyGID64();
00128   myNumIDs_ = theMap.NumMyElements();
00129   nonlocalCoefs_.resize(theNumVectors);
00130 }
00131 
00132 //----------------------------------------------------------------------------
00133 Epetra_FEVector::Epetra_FEVector(const Epetra_FEVector& source)
00134   : Epetra_MultiVector(source),
00135     myFirstID_(0),
00136     myNumIDs_(0),
00137 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00138     nonlocalIDs_int_(),
00139 #endif
00140 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00141     nonlocalIDs_LL_(),
00142 #endif
00143     nonlocalElementSize_(),
00144     nonlocalCoefs_(),
00145     nonlocalMap_(0),
00146     exporter_(0),
00147     nonlocalVector_(0),
00148     ignoreNonLocalEntries_(source.ignoreNonLocalEntries_)
00149 {
00150   *this = source;
00151 }
00152 
00153 //----------------------------------------------------------------------------
00154 Epetra_FEVector::~Epetra_FEVector()
00155 {
00156   destroyNonlocalData();
00157   destroyNonlocalMapAndExporter();
00158 
00159   nonlocalCoefs_.clear();
00160 }
00161 
00162 //----------------------------------------------------------------------------
00163 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00164 int Epetra_FEVector::SumIntoGlobalValues(int numIDs, const int* GIDs,
00165                        const double* values,
00166                                          int vectorIndex)
00167 {
00168   return( inputValues( numIDs, GIDs, values, true, vectorIndex) );
00169 }
00170 #endif
00171 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00172 int Epetra_FEVector::SumIntoGlobalValues(int numIDs, const long long* GIDs,
00173                        const double* values,
00174                                          int vectorIndex)
00175 {
00176   return( inputValues( numIDs, GIDs, values, true, vectorIndex) );
00177 }
00178 #endif
00179 //----------------------------------------------------------------------------
00180 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00181 int Epetra_FEVector::SumIntoGlobalValues(const Epetra_IntSerialDenseVector& GIDs,
00182                        const Epetra_SerialDenseVector& values,
00183                                          int vectorIndex)
00184 {
00185   if (GIDs.Length() != values.Length()) {
00186     return(-1);
00187   }
00188 
00189   return( inputValues( GIDs.Length(), GIDs.Values(), values.Values(), true,
00190                        vectorIndex ) );
00191 }
00192 #endif
00193 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00194 int Epetra_FEVector::SumIntoGlobalValues(const Epetra_LongLongSerialDenseVector& GIDs,
00195                        const Epetra_SerialDenseVector& values,
00196                                          int vectorIndex)
00197 {
00198   if (GIDs.Length() != values.Length()) {
00199     return(-1);
00200   }
00201 
00202   return( inputValues( GIDs.Length(), GIDs.Values(), values.Values(), true,
00203                        vectorIndex ) );
00204 }
00205 #endif
00206 //----------------------------------------------------------------------------
00207 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00208 int Epetra_FEVector::SumIntoGlobalValues(int numIDs, const int* GIDs,
00209            const int* numValuesPerID,
00210                        const double* values,
00211                                          int vectorIndex)
00212 {
00213   return( inputValues( numIDs, GIDs, numValuesPerID, values, true,
00214                        vectorIndex) );
00215 }
00216 #endif
00217 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00218 int Epetra_FEVector::SumIntoGlobalValues(int numIDs, const long long* GIDs,
00219            const int* numValuesPerID,
00220                        const double* values,
00221                                          int vectorIndex)
00222 {
00223   return( inputValues( numIDs, GIDs, numValuesPerID, values, true,
00224                        vectorIndex) );
00225 }
00226 #endif
00227 //----------------------------------------------------------------------------
00228 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00229 int Epetra_FEVector::ReplaceGlobalValues(int numIDs, const int* GIDs,
00230                        const double* values,
00231                                          int vectorIndex)
00232 {
00233   return( inputValues( numIDs, GIDs, values, false,
00234                        vectorIndex) );
00235 }
00236 #endif
00237 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00238 int Epetra_FEVector::ReplaceGlobalValues(int numIDs, const long long* GIDs,
00239                        const double* values,
00240                                          int vectorIndex)
00241 {
00242   return( inputValues( numIDs, GIDs, values, false,
00243                        vectorIndex) );
00244 }
00245 #endif
00246 //----------------------------------------------------------------------------
00247 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00248 int Epetra_FEVector::ReplaceGlobalValues(const Epetra_IntSerialDenseVector& GIDs,
00249                        const Epetra_SerialDenseVector& values,
00250                                          int vectorIndex)
00251 {
00252   if (GIDs.Length() != values.Length()) {
00253     return(-1);
00254   }
00255 
00256   return( inputValues( GIDs.Length(), GIDs.Values(), values.Values(), false,
00257                        vectorIndex) );
00258 }
00259 #endif
00260 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00261 int Epetra_FEVector::ReplaceGlobalValues(const Epetra_LongLongSerialDenseVector& GIDs,
00262                        const Epetra_SerialDenseVector& values,
00263                                          int vectorIndex)
00264 {
00265   if (GIDs.Length() != values.Length()) {
00266     return(-1);
00267   }
00268 
00269   return( inputValues( GIDs.Length(), GIDs.Values(), values.Values(), false,
00270                        vectorIndex) );
00271 }
00272 #endif
00273 //----------------------------------------------------------------------------
00274 template<typename int_type>
00275 int Epetra_FEVector::inputValues(int numIDs,
00276                                  const int_type* GIDs,
00277                                  const double* values,
00278                                  bool suminto,
00279                                  int vectorIndex)
00280 {
00281  //Important note!! This method assumes that there is only 1 point
00282  //associated with each element (GID), and writes to offset 0 in that
00283  //GID's block.
00284 
00285   for(int i=0; i<numIDs; ++i) {
00286     if (Map().MyGID(GIDs[i])) {
00287       if (suminto) {
00288         SumIntoGlobalValue(GIDs[i], 0, vectorIndex, values[i]);
00289       }
00290       else {
00291         ReplaceGlobalValue(GIDs[i], 0, vectorIndex, values[i]);
00292       }
00293     }
00294     else {
00295       if (!ignoreNonLocalEntries_) {
00296         EPETRA_CHK_ERR( inputNonlocalValue(GIDs[i], values[i], suminto,
00297               vectorIndex) );
00298       }
00299     }
00300   }
00301 
00302   return(0);
00303 }
00304 
00305 //----------------------------------------------------------------------------
00306 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00307 int Epetra_FEVector::ReplaceGlobalValues(int numIDs, const int* GIDs,
00308            const int* numValuesPerID,
00309                        const double* values,
00310                                          int vectorIndex)
00311 {
00312   return( inputValues( numIDs, GIDs, numValuesPerID, values, false,
00313                        vectorIndex) );
00314 }
00315 #endif
00316 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00317 int Epetra_FEVector::ReplaceGlobalValues(int numIDs, const long long* GIDs,
00318            const int* numValuesPerID,
00319                        const double* values,
00320                                          int vectorIndex)
00321 {
00322   return( inputValues( numIDs, GIDs, numValuesPerID, values, false,
00323                        vectorIndex) );
00324 }
00325 #endif
00326 //----------------------------------------------------------------------------
00327 template<typename int_type>
00328 int Epetra_FEVector::inputValues(int numIDs,
00329                                  const int_type* GIDs,
00330                                  const int* numValuesPerID,
00331                                  const double* values,
00332                                  bool suminto,
00333                                  int vectorIndex)
00334 {
00335   if(!Map().template GlobalIndicesIsType<int_type>())
00336   throw ReportError("Epetra_FEVector::inputValues mismatch between argument types (int/long long) and map type.", -1);
00337 
00338   int offset=0;
00339   for(int i=0; i<numIDs; ++i) {
00340     int numValues = numValuesPerID[i];
00341     if (Map().MyGID(GIDs[i])) {
00342       if (suminto) {
00343         for(int j=0; j<numValues; ++j) {
00344           SumIntoGlobalValue(GIDs[i], j, vectorIndex, values[offset+j]);
00345         }
00346       }
00347       else {
00348         for(int j=0; j<numValues; ++j) {
00349           ReplaceGlobalValue(GIDs[i], j, vectorIndex, values[offset+j]);
00350         }
00351       }
00352     }
00353     else {
00354       if (!ignoreNonLocalEntries_) {
00355         EPETRA_CHK_ERR( inputNonlocalValues(GIDs[i], numValues,
00356               &(values[offset]), suminto,
00357               vectorIndex) );
00358       }
00359     }
00360     offset += numValues;
00361   }
00362 
00363   return(0);
00364 }
00365 
00366 //----------------------------------------------------------------------------
00367 template<typename int_type>
00368 int Epetra_FEVector::inputNonlocalValue(int_type GID, double value, bool suminto,
00369                                         int vectorIndex)
00370 {
00371   return inputNonlocalValues(GID, 1, &value, suminto, vectorIndex);
00372 }
00373 
00374 //----------------------------------------------------------------------------
00375 template<typename int_type>
00376 int Epetra_FEVector::inputNonlocalValues(int_type GID, int numValues,
00377                                          const double* values, bool suminto,
00378                                          int vectorIndex)
00379 {
00380   if(!Map().template GlobalIndicesIsType<int_type>())
00381   throw ReportError("Epetra_FEVector::inputValues mismatch between argument types (int/long long) and map type.", -1);
00382 
00383 
00384   //find offset of GID in nonlocalIDs_var
00385 
00386   std::vector<int_type>& nonlocalIDs_var = nonlocalIDs<int_type>();
00387 
00388   typename std::vector<int_type>::iterator it = std::lower_bound(nonlocalIDs_var.begin(), nonlocalIDs_var.end(), GID);
00389   int offset = (int) (it - nonlocalIDs_var.begin());
00390   int insertPoint = offset;
00391   if (it == nonlocalIDs_var.end() || *it != GID) {
00392     offset = -1;
00393   }
00394 
00395   int elemSize = Map().MaxElementSize();
00396   if (offset >= 0) {
00397     //if offset >= 0 (meaning GID was found)
00398     //  put value in nonlocalCoefs_[vectorIndex][offset*elemSize]
00399 
00400     if (numValues != nonlocalElementSize_[offset]) {
00401       std::cerr << "Epetra_FEVector ERROR: block-size for GID " << GID << " is "
00402      << numValues<<" which doesn't match previously set block-size of "
00403      << nonlocalElementSize_[offset] << std::endl;
00404       return(-1);
00405     }
00406 
00407     offset = offset*elemSize;
00408 
00409     if (suminto) {
00410       for(int j=0; j<numValues; ++j) {
00411         nonlocalCoefs_[vectorIndex][offset+j] += values[j];
00412       }
00413     }
00414     else {
00415       for(int j=0; j<numValues; ++j) {
00416         nonlocalCoefs_[vectorIndex][offset+j] = values[j];
00417       }
00418     }
00419   }
00420   else {
00421     //else
00422     //  insert GID in nonlocalIDs_
00423     //  insert numValues   in nonlocalElementSize_
00424     //  insert values in nonlocalCoefs_
00425 
00426     nonlocalIDs_var.insert(it, GID);
00427     nonlocalElementSize_.insert(nonlocalElementSize_.begin()+insertPoint, numValues);
00428 
00429     //to keep nonlocalCoefs_[i] the same length for each vector in the multi-
00430     //vector, we'll insert positions for each vector even though values are
00431     //only being set for one of them...
00432     for(int i=0; i<NumVectors(); ++i) {
00433       for(int ii=0; ii<elemSize; ++ii) {
00434         nonlocalCoefs_[i].insert(nonlocalCoefs_[i].begin()+insertPoint*elemSize+ii, 0.0);
00435       }
00436     }
00437 
00438     for(int j=0; j<numValues; ++j) {
00439       nonlocalCoefs_[vectorIndex][insertPoint*elemSize+j] = values[j];
00440     }
00441   }
00442 
00443   return(0);
00444 }
00445 
00446 //----------------------------------------------------------------------------
00447 template<typename int_type>
00448 int Epetra_FEVector::GlobalAssemble(Epetra_CombineMode mode,
00449                                     bool reuse_map_and_exporter)
00450 {
00451   //In this method we need to gather all the non-local (overlapping) data
00452   //that's been input on each processor, into the (probably) non-overlapping
00453   //distribution defined by the map that 'this' vector was constructed with.
00454 
00455   //We don't need to do anything if there's only one processor or if
00456   //ignoreNonLocalEntries_ is true.
00457   if (Map().Comm().NumProc() < 2 || ignoreNonLocalEntries_) {
00458     return(0);
00459   }
00460 
00461   if (nonlocalMap_ == 0 || !reuse_map_and_exporter) {
00462     createNonlocalMapAndExporter<int_type>();
00463   }
00464 
00465   Epetra_MultiVector& nonlocalVector = *nonlocalVector_;
00466   nonlocalVector.PutScalar(0.0);
00467 
00468   int elemSize = Map().MaxElementSize();
00469   for(int vi=0; vi<NumVectors(); ++vi) {
00470     for(size_t i=0; i<nonlocalIDs<int_type>().size(); ++i) {
00471       for(int j=0; j<nonlocalElementSize_[i]; ++j) {
00472         nonlocalVector.ReplaceGlobalValue(nonlocalIDs<int_type>()[i], j, vi,
00473                                           nonlocalCoefs_[vi][i*elemSize+j]);
00474       }
00475     }
00476   }
00477 
00478   EPETRA_CHK_ERR( Export(nonlocalVector, *exporter_, mode) );
00479 
00480   if (reuse_map_and_exporter) {
00481     zeroNonlocalData<int_type>();
00482   }
00483   else {
00484     destroyNonlocalData();
00485   }
00486 
00487   return(0);
00488 }
00489 
00490 int Epetra_FEVector::GlobalAssemble(Epetra_CombineMode mode,
00491                                     bool reuse_map_and_exporter)
00492 {
00493   if(Map().GlobalIndicesInt())
00494 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00495     return GlobalAssemble<int>(mode, reuse_map_and_exporter);
00496 #else
00497     throw ReportError("Epetra_FEVector::GlobalAssemble: ERROR, GlobalIndicesInt but no API for it.",-1);
00498 #endif
00499 
00500   if(Map().GlobalIndicesLongLong())
00501 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00502     return GlobalAssemble<long long>(mode, reuse_map_and_exporter);
00503 #else
00504     throw ReportError("Epetra_FEVector::GlobalAssemble: ERROR, GlobalIndicesLongLong but no API for it.",-1);
00505 #endif
00506 
00507   throw ReportError("Epetra_FEVector::GlobalAssemble: Internal error, unable to determine global index type of maps", -1);
00508 }
00509 //----------------------------------------------------------------------------
00510 template<typename int_type>
00511 void Epetra_FEVector::createNonlocalMapAndExporter()
00512 {
00513   std::vector<int_type>& nonlocalIDs_var = nonlocalIDs<int_type>();
00514   delete nonlocalMap_;
00515   int_type* nlIDptr = nonlocalIDs_var.size()>0 ? &nonlocalIDs_var[0] : NULL;
00516   int* nlElSzptr = nonlocalElementSize_.size()>0 ? &nonlocalElementSize_[0] : NULL;
00517   nonlocalMap_ = new Epetra_BlockMap ((int_type) -1, (int) nonlocalIDs_var.size(), nlIDptr,
00518                                       nlElSzptr, (int_type) Map().IndexBase64(), Map().Comm());
00519   delete exporter_;
00520   exporter_ = new Epetra_Export (*nonlocalMap_, Map());
00521 
00522   delete nonlocalVector_;
00523   nonlocalVector_ = new Epetra_MultiVector (*nonlocalMap_, NumVectors());
00524 }
00525 
00526 //----------------------------------------------------------------------------
00527 void Epetra_FEVector::destroyNonlocalMapAndExporter()
00528 {
00529   delete nonlocalMap_;    nonlocalMap_ = 0;
00530   delete exporter_;       exporter_ = 0;
00531   delete nonlocalVector_; nonlocalVector_ = 0;
00532 }
00533 
00534 //----------------------------------------------------------------------------
00535 Epetra_FEVector& Epetra_FEVector::operator=(const Epetra_FEVector& source)
00536 {
00537   if (this == &source) {
00538     // Don't allow self-assignment, since the allocations and
00539     // deallocations in the code below assume that source is a
00540     // different object than *this.
00541     return *this;
00542   }
00543   // This redundantly checks for self-assignment, but the check is
00544   // inexpensive (just a pointer comparison).
00545   Epetra_MultiVector::Assign(source);
00546 
00547 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00548   nonlocalIDs_int_ = source.nonlocalIDs_int_;
00549 #endif
00550 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00551   nonlocalIDs_LL_ = source.nonlocalIDs_LL_;
00552 #endif
00553   nonlocalElementSize_ = source.nonlocalElementSize_;
00554   nonlocalCoefs_ = source.nonlocalCoefs_;
00555 
00556   return(*this);
00557 }
00558 
00559 //----------------------------------------------------------------------------
00560 template<typename int_type>
00561 void Epetra_FEVector::zeroNonlocalData()
00562 {
00563   if (nonlocalIDs<int_type>().size() > 0) {
00564     int maxelemSize = Map().MaxElementSize();
00565     for(int vi=0; vi<NumVectors(); ++vi) {
00566       for(size_t i=0; i<nonlocalIDs<int_type>().size(); ++i) {
00567         int elemSize = nonlocalElementSize_[i];
00568         for(int j=0; j<elemSize; ++j) {
00569           nonlocalCoefs_[vi][i*maxelemSize+j] = 0.0;
00570         }
00571       }
00572     }
00573   }
00574 }
00575 
00576 //----------------------------------------------------------------------------
00577 void Epetra_FEVector::destroyNonlocalData()
00578 {
00579 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00580   nonlocalIDs_int_.clear();
00581 #endif
00582 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00583   nonlocalIDs_LL_.clear();
00584 #endif
00585   nonlocalElementSize_.clear();
00586 
00587   if (nonlocalCoefs_.size() > 0) {
00588     for(int i=0; i<NumVectors(); ++i) {
00589       nonlocalCoefs_[i].clear();
00590     }
00591   }
00592 }
00593 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines