Epetra Package Browser (Single Doxygen Collection) Development
Epetra_IntVector.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_IntVector.h"
00045 #include "Epetra_Map.h"
00046 #include "Epetra_Comm.h"
00047 
00048 //=============================================================================
00049 Epetra_IntVector::Epetra_IntVector(const Epetra_BlockMap& map, bool zeroOut)
00050   : Epetra_DistObject(map, "Epetra::IntVector"),
00051     Values_(0),
00052     UserAllocated_(false),
00053     Allocated_(false)
00054 {
00055   AllocateForCopy();
00056   if(zeroOut) PutValue(0); // Zero out values
00057 }
00058 //=============================================================================
00059 Epetra_IntVector::Epetra_IntVector(const Epetra_IntVector& Source)
00060   : Epetra_DistObject(Source),
00061     Values_(0),
00062     UserAllocated_(false),
00063     Allocated_(false)
00064 {
00065   AllocateForCopy();
00066   DoCopy(Source.Values_);
00067 }
00068 //=============================================================================
00069 Epetra_IntVector::Epetra_IntVector(Epetra_DataAccess CV, const Epetra_BlockMap& map, int *V)
00070   : Epetra_DistObject(map, "Epetra::IntVector"),
00071     Values_(0),
00072     UserAllocated_(false),
00073     Allocated_(false)
00074 {
00075   if (CV==Copy) {
00076     AllocateForCopy();
00077     DoCopy(V);
00078   }
00079   else {
00080     AllocateForView();
00081     DoView(V);
00082   }
00083 }
00084 //=========================================================================
00085 Epetra_IntVector::~Epetra_IntVector(){
00086 
00087 
00088  if (Allocated_ && (!UserAllocated_) && Values_!=0) delete [] Values_;
00089 }
00090 
00091 //=========================================================================
00092 int Epetra_IntVector::AllocateForCopy()
00093 {
00094 
00095   if (Allocated_) return(0);
00096 
00097   int myLength = MyLength();
00098   if (myLength>0)
00099     Values_ = new int[myLength];
00100   else
00101     Values_ = 0;
00102 
00103   Allocated_ = true;
00104   UserAllocated_ = false;
00105   return(0);
00106 }
00107 
00108 //=========================================================================
00109 int Epetra_IntVector::DoCopy(int * V)
00110 {
00111   int iend = MyLength();
00112   for (int i=0; i<iend; i++) Values_[i] = V[i];
00113 
00114   return(0);
00115 }
00116 //=========================================================================
00117 int Epetra_IntVector::AllocateForView()
00118 {
00119 
00120   Allocated_ = true;
00121   UserAllocated_ = true;
00122 
00123   return(0);
00124 }
00125 
00126 //=========================================================================
00127 int Epetra_IntVector::DoView(int * V)
00128 {
00129 
00130   Values_ = V;
00131 
00132   return(0);
00133 }
00134 //=============================================================================
00135 int Epetra_IntVector::ExtractCopy(int *V) const {
00136 
00137   int iend = MyLength();
00138   for (int i=0; i<iend; i++) V[i] = Values_[i];
00139   return(0);
00140 }
00141 
00142 //=============================================================================
00143 int Epetra_IntVector::ExtractView(int **V) const
00144 {
00145   *V = Values_;
00146   return(0);
00147 }
00148 
00149 //=============================================================================
00150 int Epetra_IntVector::PutValue(int Value) {
00151   int iend = MyLength();
00152   for (int i=0; i<iend; i++) Values_[i] = Value;
00153   return(0);
00154 }
00155 //=============================================================================
00156 int Epetra_IntVector::MaxValue() {
00157 
00158   int result = -2000000000; // Negative 2 billion is close to smallest 32 bit int
00159   int iend = MyLength();
00160   if (iend>0) result = Values_[0];
00161   for (int i=0; i<iend; i++) result = EPETRA_MAX(result, Values_[i]);
00162   int globalResult;
00163   this->Comm().MaxAll(&result, &globalResult, 1);
00164   return(globalResult);
00165 }
00166 //=============================================================================
00167 int Epetra_IntVector::MinValue() {
00168 
00169   int result = 2000000000; // 2 billion is close to largest 32 bit int
00170   int iend = MyLength();
00171   if (iend>0) result = Values_[0];
00172   for (int i=0; i<iend; i++) result = EPETRA_MIN(result, Values_[i]);
00173   int globalResult;
00174   this->Comm().MinAll(&result, &globalResult, 1);
00175   return(globalResult);
00176 }
00177 //========================================================================
00178 Epetra_IntVector& Epetra_IntVector::operator = (const Epetra_IntVector& V) {
00179 
00180 
00181   if (MyLength() != V.MyLength())
00182     throw ReportError("Length of IntVectors incompatible in Assign.  The this IntVector has MyLength = " + toString(MyLength())
00183           + ".  The V IntVector has MyLength = " + toString(V.MyLength()), -1);
00184 
00185   int iend = MyLength();
00186   for (int i=0; i<iend; i++) Values_[i] =V[i];
00187   return(*this);
00188 }
00189 
00190 void Epetra_IntVector::Print(std::ostream& os) const {
00191   int MyPID = Map().Comm().MyPID();
00192   int NumProc = Map().Comm().NumProc();
00193 
00194   for (int iproc=0; iproc < NumProc; iproc++) {
00195     if (MyPID==iproc) {
00196       int NumMyElements1 =Map(). NumMyElements();
00197       int MaxElementSize1 = Map().MaxElementSize();
00198       int * MyGlobalElements1_int = 0;
00199       long long * MyGlobalElements1_LL = 0;
00200       if(Map().GlobalIndicesInt()) {
00201 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00202         MyGlobalElements1_int = Map().MyGlobalElements();
00203 #else
00204         throw ReportError("Epetra_IntVector::Print: Global indices int but no API for it.",-1);
00205 #endif
00206       }
00207       else if(Map().GlobalIndicesLongLong()) {
00208 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00209         MyGlobalElements1_LL = Map().MyGlobalElements64();
00210 #else
00211         throw ReportError("Epetra_IntVector::Print: Global indices long long but no API for it.",-1);
00212 #endif
00213       }
00214       int * FirstPointInElementList1=0;
00215       if (MaxElementSize1!=1) FirstPointInElementList1 = Map().FirstPointInElementList();
00216 
00217       if (MyPID==0) {
00218   os.width(8);
00219   os <<  "     MyPID"; os << "    ";
00220   os.width(12);
00221   if (MaxElementSize1==1)
00222     os <<  "GID  ";
00223   else
00224     os <<  "     GID/Point";
00225   os.width(20);
00226   os <<  "Value  ";
00227   os << std::endl;
00228       }
00229       for (int i=0; i < NumMyElements1; i++) {
00230   for (int ii=0; ii< Map().ElementSize(i); ii++) {
00231     int iii;
00232     os.width(10);
00233     os <<  MyPID; os << "    ";
00234     os.width(10);
00235     if (MaxElementSize1==1) {
00236         if(MyGlobalElements1_int)
00237         os << MyGlobalElements1_int[i] << "    ";
00238         if(MyGlobalElements1_LL)
00239         os << MyGlobalElements1_LL[i] << "    ";
00240       iii = i;
00241     }
00242     else {
00243         if(MyGlobalElements1_int)
00244         os <<  MyGlobalElements1_int[i]<< "/" << ii << "    ";
00245         if(MyGlobalElements1_LL)
00246         os <<  MyGlobalElements1_LL[i]<< "/" << ii << "    ";
00247       iii = FirstPointInElementList1[i]+ii;
00248     }
00249         os.width(20);
00250         os <<  Values_[iii];
00251     os << std::endl;
00252   }
00253       }
00254       os << std::flush;
00255     }
00256 
00257     // Do a few global ops to give I/O a chance to complete
00258     Map().Comm().Barrier();
00259     Map().Comm().Barrier();
00260     Map().Comm().Barrier();
00261   }
00262   return;
00263 }
00264 //=========================================================================
00265 int Epetra_IntVector::CheckSizes(const Epetra_SrcDistObject& Source)
00266 {
00267   (void)Source;
00268   return(0);
00269 }
00270 
00271 //=========================================================================
00272 int Epetra_IntVector::CopyAndPermute(const Epetra_SrcDistObject& Source,
00273                                      int NumSameIDs,
00274                                      int NumPermuteIDs,
00275                                      int * PermuteToLIDs,
00276                                      int *PermuteFromLIDs,
00277                                      const Epetra_OffsetIndex * Indexor)
00278 {
00279   (void)Indexor;
00280   const Epetra_IntVector & A = dynamic_cast<const Epetra_IntVector &>(Source);
00281 
00282   int * From;
00283   A.ExtractView(&From);
00284   int *To = Values_;
00285 
00286   int * ToFirstPointInElementList = 0;
00287   int * FromFirstPointInElementList = 0;
00288   int * FromElementSizeList = 0;
00289   int MaxElementSize = Map().MaxElementSize();
00290   bool ConstantElementSize = Map().ConstantElementSize();
00291 
00292   if (!ConstantElementSize) {
00293     ToFirstPointInElementList =   Map().FirstPointInElementList();
00294     FromFirstPointInElementList = A.Map().FirstPointInElementList();
00295     FromElementSizeList = A.Map().ElementSizeList();
00296   }
00297   int j, jj, jjj, k;
00298 
00299   int NumSameEntries;
00300 
00301   bool Case1 = false;
00302   bool Case2 = false;
00303   // bool Case3 = false;
00304 
00305   if (MaxElementSize==1) {
00306     Case1 = true;
00307     NumSameEntries = NumSameIDs;
00308   }
00309   else if (ConstantElementSize) {
00310     Case2 = true;
00311     NumSameEntries = NumSameIDs * MaxElementSize;
00312   }
00313   else {
00314     // Case3 = true;
00315     NumSameEntries = FromFirstPointInElementList[NumSameIDs];
00316   }
00317 
00318   // Short circuit for the case where the source and target vector is the same.
00319   if (To==From) NumSameEntries = 0;
00320 
00321   // Do copy first
00322   if (NumSameIDs>0)
00323     if (To!=From) {
00324   for (j=0; j<NumSameEntries; j++)
00325     To[j] = From[j];
00326     }
00327   // Do local permutation next
00328   if (NumPermuteIDs>0) {
00329 
00330     // Point entry case
00331     if (Case1) {
00332 
00333       for (j=0; j<NumPermuteIDs; j++)
00334   To[PermuteToLIDs[j]] = From[PermuteFromLIDs[j]];
00335     }
00336     // constant element size case
00337     else if (Case2) {
00338 
00339       for (j=0; j<NumPermuteIDs; j++) {
00340   jj = MaxElementSize*PermuteToLIDs[j];
00341   jjj = MaxElementSize*PermuteFromLIDs[j];
00342     for (k=0; k<MaxElementSize; k++)
00343       To[jj+k] = From[jjj+k];
00344       }
00345     }
00346 
00347     // variable element size case
00348     else {
00349 
00350       for (j=0; j<NumPermuteIDs; j++) {
00351   jj = ToFirstPointInElementList[PermuteToLIDs[j]];
00352   jjj = FromFirstPointInElementList[PermuteFromLIDs[j]];
00353   int ElementSize = FromElementSizeList[PermuteFromLIDs[j]];
00354     for (k=0; k<ElementSize; k++)
00355       To[jj+k] = From[jjj+k];
00356       }
00357     }
00358   }
00359   return(0);
00360 }
00361 
00362 //=========================================================================
00363 int Epetra_IntVector::PackAndPrepare(const Epetra_SrcDistObject & Source,
00364                                      int NumExportIDs,
00365                                      int * ExportLIDs,
00366                                      int & LenExports,
00367                                      char * & Exports,
00368                                      int & SizeOfPacket,
00369                                      int * Sizes,
00370                                      bool & VarSizes,
00371                                      Epetra_Distributor & Distor)
00372 {
00373   (void)Sizes;
00374   (void)VarSizes;
00375   (void)Distor;
00376   const Epetra_IntVector & A = dynamic_cast<const Epetra_IntVector &>(Source);
00377 
00378   int j, jj, k;
00379 
00380   int  * From;
00381   A.ExtractView(&From);
00382   int MaxElementSize = Map().MaxElementSize();
00383   bool ConstantElementSize = Map().ConstantElementSize();
00384 
00385   int * FromFirstPointInElementList = 0;
00386   int * FromElementSizeList = 0;
00387 
00388   if (!ConstantElementSize) {
00389     FromFirstPointInElementList = A.Map().FirstPointInElementList();
00390     FromElementSizeList = A.Map().ElementSizeList();
00391   }
00392 
00393   SizeOfPacket = MaxElementSize * (int)sizeof(int);
00394 
00395   if(NumExportIDs*SizeOfPacket>LenExports) {
00396     if (LenExports>0) delete [] Exports;
00397     LenExports = NumExportIDs*SizeOfPacket;
00398     Exports = new char[LenExports];
00399   }
00400 
00401   int * ptr;
00402 
00403   if (NumExportIDs>0) {
00404     ptr = (int *) Exports;
00405 
00406     // Point entry case
00407     if (MaxElementSize==1) for (j=0; j<NumExportIDs; j++) *ptr++ = From[ExportLIDs[j]];
00408 
00409     // constant element size case
00410     else if (ConstantElementSize) {
00411 
00412       for (j=0; j<NumExportIDs; j++) {
00413   jj = MaxElementSize*ExportLIDs[j];
00414     for (k=0; k<MaxElementSize; k++)
00415       *ptr++ = From[jj+k];
00416       }
00417     }
00418 
00419     // variable element size case
00420     else {
00421 
00422       int thisSizeOfPacket = MaxElementSize;
00423       for (j=0; j<NumExportIDs; j++) {
00424   ptr = (int *) Exports + j*thisSizeOfPacket;
00425   jj = FromFirstPointInElementList[ExportLIDs[j]];
00426   int ElementSize = FromElementSizeList[ExportLIDs[j]];
00427     for (k=0; k<ElementSize; k++)
00428       *ptr++ = From[jj+k];
00429       }
00430     }
00431   }
00432 
00433   return(0);
00434 }
00435 
00436 //=========================================================================
00437 int Epetra_IntVector::UnpackAndCombine(const Epetra_SrcDistObject & Source,
00438                                        int NumImportIDs,
00439                                        int * ImportLIDs,
00440                                        int LenImports,
00441                                        char * Imports,
00442                                        int & SizeOfPacket,
00443                                        Epetra_Distributor & Distor,
00444                                        Epetra_CombineMode CombineMode,
00445                                        const Epetra_OffsetIndex * Indexor)
00446 {
00447   (void)Source;
00448   (void)LenImports;
00449   (void)SizeOfPacket;
00450   (void)Distor;
00451   (void)Indexor;
00452   int j, jj, k;
00453 
00454   if(    CombineMode != Add
00455       && CombineMode != Zero
00456       && CombineMode != Insert
00457       && CombineMode != Average
00458       && CombineMode != AbsMax )
00459     EPETRA_CHK_ERR(-1); //Unsupported CombinedMode, will default to Zero
00460 
00461   if (NumImportIDs<=0) return(0);
00462 
00463   int * To = Values_;
00464   int MaxElementSize = Map().MaxElementSize();
00465   bool ConstantElementSize = Map().ConstantElementSize();
00466 
00467   int * ToFirstPointInElementList = 0;
00468   int * ToElementSizeList = 0;
00469 
00470   if (!ConstantElementSize) {
00471     ToFirstPointInElementList = Map().FirstPointInElementList();
00472     ToElementSizeList = Map().ElementSizeList();
00473   }
00474 
00475   int * ptr;
00476   // Unpack it...
00477 
00478   ptr = (int *) Imports;
00479 
00480   // Point entry case
00481   if (MaxElementSize==1) {
00482 
00483       if (CombineMode==Add)
00484   for (j=0; j<NumImportIDs; j++) To[ImportLIDs[j]] += *ptr++; // Add to existing value
00485       else if(CombineMode==Insert)
00486   for (j=0; j<NumImportIDs; j++) To[ImportLIDs[j]] = *ptr++;
00487       else if(CombineMode==AbsMax)
00488         for (j=0; j<NumImportIDs; j++) {
00489     To[ImportLIDs[j]] = EPETRA_MAX( To[ImportLIDs[j]],std::abs(*ptr));
00490     ptr++;
00491   }
00492       // Note:  The following form of averaging is not a true average if more that one value is combined.
00493       //        This might be an issue in the future, but we leave this way for now.
00494       else if(CombineMode==Average)
00495   for (j=0; j<NumImportIDs; j++) {To[ImportLIDs[j]] += *ptr++; To[ImportLIDs[j]] /= 2;}
00496   }
00497 
00498   // constant element size case
00499 
00500   else if (ConstantElementSize) {
00501 
00502     if (CombineMode==Add) {
00503       for (j=0; j<NumImportIDs; j++) {
00504   jj = MaxElementSize*ImportLIDs[j];
00505     for (k=0; k<MaxElementSize; k++)
00506       To[jj+k] += *ptr++; // Add to existing value
00507       }
00508     }
00509     else if(CombineMode==Insert) {
00510       for (j=0; j<NumImportIDs; j++) {
00511   jj = MaxElementSize*ImportLIDs[j];
00512     for (k=0; k<MaxElementSize; k++)
00513       To[jj+k] = *ptr++;
00514       }
00515     }
00516     else if(CombineMode==AbsMax) {
00517       for (j=0; j<NumImportIDs; j++) {
00518   jj = MaxElementSize*ImportLIDs[j];
00519   for (k=0; k<MaxElementSize; k++) {
00520       To[jj+k] = EPETRA_MAX( To[jj+k], std::abs(*ptr));
00521       ptr++;
00522   }
00523       }
00524     }
00525     // Note:  The following form of averaging is not a true average if more that one value is combined.
00526     //        This might be an issue in the future, but we leave this way for now.
00527     else if(CombineMode==Average) {
00528       for (j=0; j<NumImportIDs; j++) {
00529   jj = MaxElementSize*ImportLIDs[j];
00530     for (k=0; k<MaxElementSize; k++)
00531       { To[jj+k] += *ptr++; To[jj+k] /= 2;}
00532       }
00533     }
00534   }
00535 
00536   // variable element size case
00537 
00538   else {
00539 
00540     int thisSizeOfPacket = MaxElementSize;
00541 
00542     if (CombineMode==Add) {
00543       for (j=0; j<NumImportIDs; j++) {
00544   ptr = (int *) Imports + j*thisSizeOfPacket;
00545   jj = ToFirstPointInElementList[ImportLIDs[j]];
00546   int ElementSize = ToElementSizeList[ImportLIDs[j]];
00547     for (k=0; k<ElementSize; k++)
00548       To[jj+k] += *ptr++; // Add to existing value
00549       }
00550     }
00551     else  if(CombineMode==Insert){
00552       for (j=0; j<NumImportIDs; j++) {
00553   ptr = (int *) Imports + j*thisSizeOfPacket;
00554   jj = ToFirstPointInElementList[ImportLIDs[j]];
00555   int ElementSize = ToElementSizeList[ImportLIDs[j]];
00556     for (k=0; k<ElementSize; k++)
00557       To[jj+k] = *ptr++;
00558       }
00559     }
00560     else  if(CombineMode==AbsMax){
00561       for (j=0; j<NumImportIDs; j++) {
00562   ptr = (int *) Imports + j*thisSizeOfPacket;
00563   jj = ToFirstPointInElementList[ImportLIDs[j]];
00564   int ElementSize = ToElementSizeList[ImportLIDs[j]];
00565   for (k=0; k<ElementSize; k++) {
00566       To[jj+k] = EPETRA_MAX( To[jj+k], std::abs(*ptr));
00567       ptr++;
00568   }
00569       }
00570     }
00571     // Note:  The following form of averaging is not a true average if more that one value is combined.
00572     //        This might be an issue in the future, but we leave this way for now.
00573     else if(CombineMode==Average) {
00574       for (j=0; j<NumImportIDs; j++) {
00575   ptr = (int *) Imports + j*thisSizeOfPacket;
00576   jj = ToFirstPointInElementList[ImportLIDs[j]];
00577   int ElementSize = ToElementSizeList[ImportLIDs[j]];
00578     for (k=0; k<ElementSize; k++)
00579       { To[jj+k] += *ptr++; To[jj+k] /= 2;}
00580       }
00581     }
00582   }
00583 
00584   return(0);
00585 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines