Epetra Package Browser (Single Doxygen Collection) Development
Epetra_BasicRowMatrix.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_BasicRowMatrix.h"
00044 #include "Epetra_Map.h"
00045 #include "Epetra_Import.h"
00046 #include "Epetra_Export.h"
00047 #include "Epetra_Vector.h"
00048 #include "Epetra_MultiVector.h"
00049 #include "Epetra_Comm.h"
00050 #include "Epetra_Util.h"
00051 #include "Epetra_IntSerialDenseVector.h"
00052 #include "Epetra_RowMatrix.h"
00053 #include "Epetra_CrsMatrix.h"
00054 
00055 
00056 //==============================================================================
00057 Epetra_BasicRowMatrix::Epetra_BasicRowMatrix(const Epetra_Comm & comm) 
00058   : Comm_(comm.Clone()),
00059     OperatorDomainMap_(Epetra_Map(0,0,comm)),
00060     OperatorRangeMap_(Epetra_Map(0,0,comm)),
00061     RowMatrixRowMap_(Epetra_Map(0,0,comm)),
00062     RowMatrixColMap_(Epetra_Map(0,0,comm)),
00063     NumMyNonzeros_(0),
00064     NumGlobalNonzeros_(0),
00065     MaxNumEntries_(0),
00066     NormInf_(0.0),
00067     NormOne_(0.0),
00068     UseTranspose_(false),
00069     HasNormInf_(true),
00070     LowerTriangular_(true),
00071     UpperTriangular_(true),
00072     HaveStructureConstants_(false),
00073     HaveNumericConstants_(false),
00074     HaveMaps_(false),
00075     ImportVector_(0),
00076     ExportVector_(0),
00077     Importer_(0),
00078     Exporter_(0)
00079 {
00080   SetLabel("Epetra::BasicRowMatrix");
00081 }
00082 
00083 //==============================================================================
00084 Epetra_BasicRowMatrix::~Epetra_BasicRowMatrix(){
00085 
00086   if (ImportVector_!=0) delete ImportVector_;
00087   ImportVector_=0;
00088   if (ExportVector_!=0) delete ExportVector_;
00089   ExportVector_=0;
00090   if (Importer_!=0) delete Importer_;
00091   Importer_=0;
00092   if (Exporter_!=0) delete Exporter_;
00093   Exporter_=0;
00094   delete Comm_;
00095 }
00096 
00097 //==============================================================================
00098 void Epetra_BasicRowMatrix::SetMaps(const Epetra_Map & RowMap, const Epetra_Map & ColMap) {
00099 
00100   SetMaps(RowMap, ColMap, RowMap, RowMap);
00101 }
00102 
00103 //==============================================================================
00104 void Epetra_BasicRowMatrix::SetMaps(const Epetra_Map & RowMap, const Epetra_Map & ColMap, 
00105       const Epetra_Map & DomainMap, const Epetra_Map & RangeMap) {
00106 
00107   RowMatrixRowMap_ = RowMap;
00108   RowMatrixColMap_ = ColMap;
00109   OperatorDomainMap_ = DomainMap;
00110   OperatorRangeMap_ = RangeMap;
00111   HaveMaps_ = true;
00112   HaveStructureConstants_ = false;
00113   HaveNumericConstants_ = false;
00114 
00115   if (!OperatorDomainMap().UniqueGIDs()) throw RowMatrixRowMap().ReportError("At least one GID is repeated in domain map. Domain and range maps must have unique GIDs", -1);
00116   if (!OperatorRangeMap().UniqueGIDs()) throw RowMatrixRowMap().ReportError("At least one GID is repeated in range map. Domain and range maps must have unique GIDs", -2);
00117   SetImportExport();
00118 }
00119 
00120 //==============================================================================
00121 void Epetra_BasicRowMatrix::SetImportExport() {
00122 
00123   // Check if non-trivial import/export operators
00124   if (!(RowMatrixRowMap().SameAs(OperatorRangeMap()))) 
00125     Exporter_ = new Epetra_Export(RowMatrixRowMap(), OperatorRangeMap());
00126   
00127   if (!(RowMatrixColMap().SameAs(OperatorDomainMap())))
00128     Importer_ = new Epetra_Import(RowMatrixColMap(), OperatorDomainMap());
00129 
00130   NumMyRows_ = RowMatrixRowMap().NumMyPoints();
00131   NumMyCols_ = RowMatrixColMap().NumMyPoints();
00132 }
00133 
00134 //==============================================================================
00135 void Epetra_BasicRowMatrix::ComputeStructureConstants() const {
00136   MaxNumEntries_ = 0;
00137   NumMyNonzeros_ = 0;
00138   NumGlobalNonzeros_ = 0;
00139   int NumEntries = 0;
00140   for (int i=0; i<NumMyRows_; i++) {
00141     NumMyRowEntries(i, NumEntries);
00142     NumMyNonzeros_ += NumEntries;
00143     if (NumEntries>MaxNumEntries_) MaxNumEntries_ = NumEntries;
00144   }
00145 
00146   RowMatrixRowMap().Comm().SumAll(&NumMyNonzeros_, &NumGlobalNonzeros_, 1);
00147   HaveStructureConstants_ = true;
00148 }
00149 //=============================================================================
00150 void Epetra_BasicRowMatrix::ComputeNumericConstants() const {
00151   Epetra_SerialDenseVector Values(MaxNumEntries());
00152   Epetra_IntSerialDenseVector Indices(MaxNumEntries());
00153   int NumEntries;
00154   Epetra_Vector x1(RowMatrixRowMap()); // Need temp vector for row sums
00155   Epetra_Vector x2(RowMatrixColMap()); // Need temp vector for column sums
00156   for(int i = 0; i < NumMyRows_; i++) {
00157     ExtractMyRowCopy(i, MaxNumEntries(), NumEntries, Values.Values(), Indices.Values());
00158     for(int j = 0; j < NumEntries; j++) {
00159       x1[i] += std::abs(Values[j]);
00160       x2[Indices[j]] += std::abs(Values[j]);
00161       if (Indices[j]<i) UpperTriangular_ = false;
00162       if (Indices[j]>i) LowerTriangular_ = false;
00163     }
00164   }
00165 
00166   // If we have a non-trivial exporter, we must export elements that are permuted or belong to other processors
00167   if(Exporter() != 0) {
00168     Epetra_Vector xtmp(OperatorRangeMap()); // Create temporary import vector if needed
00169     xtmp.Export(x1,*Exporter(),Add);
00170     xtmp.MaxValue(&NormInf_); // This is the NormInf
00171   }
00172   else
00173     x1.MaxValue(&NormInf_); // Find max
00174 
00175   // If we have a non-trivial importer, we must export elements that are permuted or belong to other processors
00176   if(Importer() != 0) {
00177     Epetra_Vector xtmp(OperatorDomainMap()); // Create temporary import vector if needed
00178     xtmp.Export(x2,*Importer(),Add);
00179     xtmp.MaxValue(&NormOne_); // This is the NormOne
00180   }
00181   else
00182     x2.MaxValue(&NormOne_); // Find max
00183 
00184   UpdateFlops(2*NumGlobalNonzeros());
00185   HaveNumericConstants_ = true;
00186 }
00187 //=============================================================================
00188 int Epetra_BasicRowMatrix::ExtractDiagonalCopy(Epetra_Vector & Diagonal) const {
00189 
00190   if(!RowMatrixRowMap().SameAs(Diagonal.Map())) 
00191     EPETRA_CHK_ERR(-2); // Maps must be the same
00192 
00193   // Crude implementation in terms of ExtractMyRowCopy
00194 
00195   Epetra_SerialDenseVector Values(MaxNumEntries());
00196   Epetra_IntSerialDenseVector Indices(MaxNumEntries());
00197   int NumEntries;
00198 
00199   for(int i = 0; i < NumMyRows_; i++) {
00200     EPETRA_CHK_ERR(ExtractMyRowCopy(i, MaxNumEntries(), NumEntries, Values.Values(), Indices.Values()));
00201     int ii = RowMatrixRowMap().GID(i);
00202     
00203     Diagonal[i] = 0.0;
00204     for(int j = 0; j < NumEntries; j++) {
00205       if(ii == RowMatrixColMap().GID(Indices[j])) {
00206   Diagonal[i] = Values[j];
00207   break;
00208       }
00209     }
00210   }
00211   return(0);
00212 }
00213 //=============================================================================
00214 int Epetra_BasicRowMatrix::InvRowSums(Epetra_Vector & x) const {
00215   int ierr = 0;
00216   int i, j;
00217   Epetra_SerialDenseVector Values(MaxNumEntries());
00218   Epetra_IntSerialDenseVector Indices(MaxNumEntries());
00219   int NumEntries;
00220   x.PutScalar(0.0); // Make sure we sum into a vector of zeros.
00221   double * xp = (double*)x.Values();
00222   if (OperatorRangeMap().SameAs(x.Map()) && Exporter() != 0) {
00223     Epetra_Vector x_tmp(RowMatrixRowMap());
00224     x_tmp.PutScalar(0.0);
00225     double * x_tmp_p = (double*)x_tmp.Values();
00226     for (i=0; i < NumMyRows_; i++) {
00227       EPETRA_CHK_ERR(ExtractMyRowCopy(i, MaxNumEntries(), NumEntries, Values.Values(), Indices.Values()));
00228       for (j=0; j < NumEntries; j++)  x_tmp_p[i] += std::abs(Values[j]);
00229     }
00230     EPETRA_CHK_ERR(x.Export(x_tmp, *Exporter(), Add)); //Export partial row sums to x.
00231     int myLength = x.MyLength();
00232     for (i=0; i<myLength; i++) { 
00233       if (xp[i]<Epetra_MinDouble) {
00234         if (xp[i]==0.0) ierr = 1; // Set error to 1 to signal that zero rowsum found (supercedes ierr = 2)
00235         else if (ierr!=1) ierr = 2;
00236         xp[i] = Epetra_MaxDouble;
00237       }
00238       else
00239         xp[i] = 1.0/xp[i];
00240     }
00241   }
00242   else if (RowMatrixRowMap().SameAs(x.Map())) {
00243     for (i=0; i < NumMyRows_; i++) {
00244       EPETRA_CHK_ERR(ExtractMyRowCopy(i, MaxNumEntries(), NumEntries, Values.Values(), Indices.Values()));
00245       double scale = 0.0;
00246       for (j=0; j < NumEntries; j++) scale += std::abs(Values[j]);
00247       if (scale<Epetra_MinDouble) {
00248         if (scale==0.0) ierr = 1; // Set error to 1 to signal that zero rowsum found (supercedes ierr = 2)
00249         else if (ierr!=1) ierr = 2;
00250         xp[i] = Epetra_MaxDouble;
00251       }
00252       else
00253         xp[i] = 1.0/scale;
00254     }
00255   }
00256   else { // x.Map different than both RowMatrixRowMap() and OperatorRangeMap()
00257     EPETRA_CHK_ERR(-2); // The map of x must be the RowMap or RangeMap of A.
00258   }
00259   EPETRA_CHK_ERR(ierr);  
00260   UpdateFlops(NumGlobalNonzeros());
00261   return(0);
00262 }
00263 //=============================================================================
00264 int Epetra_BasicRowMatrix::LeftScale(const Epetra_Vector & x) {
00265   double *curValue;
00266   int curRowIndex, curColIndex;
00267   if(OperatorRangeMap().SameAs(x.Map()) && Exporter() != 0) {
00268     Epetra_Vector xtmp(RowMatrixRowMap());
00269     xtmp.Import(x,*Exporter(),Insert);
00270     for (int i=0; i<NumMyNonzeros_; i++) {
00271       EPETRA_CHK_ERR(ExtractMyEntryView(i, curValue, curRowIndex, curColIndex));
00272       *curValue *= xtmp[curRowIndex];
00273     }
00274   }
00275   else if (RowMatrixRowMap().SameAs(x.Map()))
00276     for (int i=0; i<NumMyNonzeros_; i++) {
00277       EPETRA_CHK_ERR(ExtractMyEntryView(i, curValue, curRowIndex, curColIndex));
00278       *curValue *= x[curRowIndex];
00279     }
00280   else {
00281     EPETRA_CHK_ERR(-2); // The Map of x must be the RowMap or RangeMap of A.
00282   }
00283   HaveNumericConstants_ = false;
00284   UpdateFlops(NumGlobalNonzeros());
00285   return(0);
00286 }
00287 //=============================================================================
00288 int Epetra_BasicRowMatrix::InvColSums(Epetra_Vector & x) const {
00289   int ierr = 0;
00290   int i, j;
00291   Epetra_SerialDenseVector Values(MaxNumEntries());
00292   Epetra_IntSerialDenseVector Indices(MaxNumEntries());
00293   int NumEntries;
00294   int MapNumMyElements = x.Map().NumMyElements();
00295   x.PutScalar(0.0); // Make sure we sum into a vector of zeros.
00296   double* xp = (double*)x.Values();
00297   if(OperatorDomainMap().SameAs(x.Map()) && Importer() != 0) {
00298     Epetra_Vector x_tmp(RowMatrixColMap());
00299     x_tmp.PutScalar(0.0);
00300     double * x_tmp_p = (double*)x_tmp.Values();
00301     for(i = 0; i < NumMyRows_; i++) {
00302       EPETRA_CHK_ERR(ExtractMyRowCopy(i, MaxNumEntries(), NumEntries, Values.Values(), Indices.Values()));
00303       for(j = 0; j < NumEntries; j++) 
00304         x_tmp_p[Indices[j]] += std::abs(Values[j]);
00305     }
00306     EPETRA_CHK_ERR(x.Export(x_tmp, *Importer(), Add)); // Fill x with partial column sums
00307   }
00308   else if(RowMatrixColMap().SameAs(x.Map())) {
00309     for(i = 0; i < NumMyRows_; i++) {
00310       EPETRA_CHK_ERR(ExtractMyRowCopy(i, MaxNumEntries(), NumEntries, Values.Values(), Indices.Values()));
00311       for(j = 0; j < NumEntries; j++) 
00312         xp[Indices[j]] += std::abs(Values[j]);
00313     }
00314   }
00315   else { //x.Map different than both RowMatrixColMap() and OperatorDomainMap()
00316     EPETRA_CHK_ERR(-2); // x must have the same distribution as the domain of A
00317   }
00318 
00319   // Invert values, don't allow them to get too large
00320   for(i = 0; i < MapNumMyElements; i++) {
00321     double scale = xp[i];
00322     if(scale < Epetra_MinDouble) {
00323       if(scale == 0.0) 
00324   ierr = 1; // Set error to 1 to signal that zero rowsum found (supercedes ierr = 2)
00325       else if(ierr != 1) 
00326   ierr = 2;
00327       xp[i] = Epetra_MaxDouble;
00328     }
00329     else
00330       xp[i] = 1.0 / scale;
00331   }
00332 
00333   EPETRA_CHK_ERR(ierr);
00334   UpdateFlops(NumGlobalNonzeros());
00335   return(0);
00336 }
00337 //=============================================================================
00338 int Epetra_BasicRowMatrix::RightScale(const Epetra_Vector & x) {
00339   double *curValue;
00340   int curRowIndex, curColIndex;
00341   if(OperatorDomainMap().SameAs(x.Map()) && Importer() != 0) {
00342     Epetra_Vector xtmp(RowMatrixColMap());
00343     xtmp.Import(x,*Importer(),Insert);
00344     for (int i=0; i<NumMyNonzeros_; i++) {
00345       EPETRA_CHK_ERR(ExtractMyEntryView(i, curValue, curRowIndex, curColIndex));
00346       *curValue *= xtmp[curColIndex];
00347     }
00348   }
00349   else if (RowMatrixColMap().SameAs(x.Map()))
00350     for (int i=0; i<NumMyNonzeros_; i++) {
00351       EPETRA_CHK_ERR(ExtractMyEntryView(i, curValue, curRowIndex, curColIndex));
00352       *curValue *= x[curColIndex];
00353     }
00354   else {
00355     EPETRA_CHK_ERR(-2); // The Map of x must be the RowMap or RangeMap of A.
00356   }
00357   HaveNumericConstants_ = false;
00358   UpdateFlops(NumGlobalNonzeros());
00359   return(0);
00360 }
00361 //=============================================================================
00362 int Epetra_BasicRowMatrix::Multiply(bool TransA, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const {
00363   //
00364   // This function forms the product Y = A * Y or Y = A' * X
00365   //
00366 
00367   Epetra_SerialDenseVector Values(MaxNumEntries());
00368   Epetra_IntSerialDenseVector Indices(MaxNumEntries());
00369   int NumEntries;
00370 
00371   int NumVectors = X.NumVectors();
00372   if (NumVectors!=Y.NumVectors()) {
00373     EPETRA_CHK_ERR(-1); // Need same number of vectors in each MV
00374   }
00375 
00376   UpdateImportVector(NumVectors); // Make sure Import and Export Vectors are compatible
00377   UpdateExportVector(NumVectors);
00378 
00379   double ** Xp = (double**) X.Pointers();
00380   double ** Yp = (double**) Y.Pointers();
00381 
00382   if (!TransA) {
00383 
00384     // If we have a non-trivial importer, we must import elements that are permuted or are on other processors
00385     if (Importer()!=0) {
00386       EPETRA_CHK_ERR(ImportVector_->Import(X, *Importer(), Insert));
00387       Xp = (double**)ImportVector_->Pointers();
00388     }
00389 
00390     // If we have a non-trivial exporter, we must export elements that are permuted or belong to other processors
00391     if (Exporter()!=0) {
00392       Yp = (double**)ExportVector_->Pointers();
00393     }
00394 
00395     // Do actual computation
00396     for(int i = 0; i < NumMyRows_; i++) {
00397       EPETRA_CHK_ERR(ExtractMyRowCopy(i, MaxNumEntries(), NumEntries, Values.Values(), Indices.Values()));
00398       for (int k=0; k<NumVectors; k++) {
00399   double sum = 0.0;
00400   for(int j = 0; j < NumEntries; j++)
00401     sum += Values[j]*Xp[k][Indices[j]];
00402   Yp[k][i] = sum;
00403       }
00404     }
00405     
00406     if (Exporter()!=0) {
00407       Y.PutScalar(0.0);  // Make sure target is zero
00408       Y.Export(*ExportVector_, *Exporter(), Add); // Fill Y with Values from export vector
00409     }
00410     // Handle case of rangemap being a local replicated map
00411     if (!OperatorRangeMap().DistributedGlobal() && Comm().NumProc()>1) EPETRA_CHK_ERR(Y.Reduce());
00412   }
00413   else { // Transpose operation
00414     
00415 
00416     // If we have a non-trivial exporter, we must import elements that are permuted or are on other processors
00417 
00418     if (Exporter()!=0) {
00419       EPETRA_CHK_ERR(ExportVector_->Import(X, *Exporter(), Insert));
00420       Xp = (double**)ExportVector_->Pointers();
00421     }
00422 
00423     // If we have a non-trivial importer, we must export elements that are permuted or belong to other processors
00424     if (Importer()!=0) {
00425       Yp = (double**)ImportVector_->Pointers();
00426       ImportVector_->PutScalar(0.0);  // Make sure target is zero
00427     }
00428     else Y.PutScalar(0.0); // Make sure target is zero
00429 
00430     // Do actual computation
00431     for(int i = 0; i < NumMyRows_; i++) {
00432       EPETRA_CHK_ERR(ExtractMyRowCopy(i, MaxNumEntries(), NumEntries, Values.Values(), Indices.Values()));
00433       for (int k=0; k<NumVectors; k++) {
00434   double xtmp = Xp[k][i];
00435   for(int j = 0; j < NumEntries; j++)
00436     Yp[k][Indices[j]] += Values[j]*xtmp;
00437       }
00438     }
00439     
00440     if (Importer()!=0) {
00441       Y.PutScalar(0.0);  // Make sure target is zero
00442       EPETRA_CHK_ERR(Y.Export(*ImportVector_, *Importer(), Add)); // Fill Y with Values from export vector
00443     }
00444     // Handle case of rangemap being a local replicated map
00445     if (!OperatorDomainMap().DistributedGlobal() && Comm().NumProc()>1)  EPETRA_CHK_ERR(Y.Reduce());
00446   }
00447 
00448   UpdateFlops(2*NumVectors*NumGlobalNonzeros());
00449   return(0);
00450 }
00451 //=======================================================================================================
00452 void Epetra_BasicRowMatrix::UpdateImportVector(int NumVectors) const {
00453   if(Importer() != 0) {
00454     if(ImportVector_ != 0) {
00455       if(ImportVector_->NumVectors() != NumVectors) {
00456      delete ImportVector_;
00457      ImportVector_= 0;
00458       }
00459     }
00460     if(ImportVector_ == 0)
00461       ImportVector_ = new Epetra_MultiVector(Importer_->TargetMap(),NumVectors); // Create import vector if needed
00462   }
00463   return;
00464 }
00465 //=======================================================================================================
00466 void Epetra_BasicRowMatrix::UpdateExportVector(int NumVectors) const {
00467   if(Exporter() != 0) {
00468     if(ExportVector_ != 0) {
00469       if(ExportVector_->NumVectors() != NumVectors) {
00470      delete ExportVector_;
00471      ExportVector_= 0;
00472       }
00473     }
00474     if(ExportVector_ == 0)
00475       ExportVector_ = new Epetra_MultiVector(Exporter_->SourceMap(),NumVectors); // Create Export vector if needed
00476   }
00477   return;
00478 }
00479 
00480 //=======================================================================================================
00481 void Epetra_BasicRowMatrix::Print(ostream& os) const {
00482 
00483   int MyPID = RowMatrixRowMap().Comm().MyPID();
00484   int NumProc = RowMatrixRowMap().Comm().NumProc();
00485 
00486   for (int iproc=0; iproc < NumProc; iproc++) {
00487     if (MyPID==iproc) {
00488       if (MyPID==0) {
00489     os <<    "\nNumber of Global Rows         = "; os << NumGlobalRows();    os << endl;
00490     os <<    "Number of Global Cols         = "; os << NumGlobalCols();    os << endl;
00491     os <<    "Number of Global Diagonals    = "; os << NumGlobalDiagonals(); os << endl;
00492   os <<    "Number of Global Nonzeros     = "; os << NumGlobalNonzeros_; os << endl;
00493       }
00494       
00495       os <<  "\nNumber of My Rows               = "; os << NumMyRows_; os << endl;
00496       os <<    "Number of My Cols               = "; os << NumMyCols_; os << endl;
00497       os <<    "Number of My Diagonals          = "; os << NumMyDiagonals(); os << endl;
00498       os <<    "Number of My Nonzeros           = "; os << NumMyNonzeros_; os << endl;
00499       os <<    "My Maximum Num Entries          = "; os << MaxNumEntries_; os << endl; os << endl;
00500       os << flush;
00501       
00502     }
00503     // Do a few global ops to give I/O a chance to complete
00504     Comm().Barrier();
00505     Comm().Barrier();
00506     Comm().Barrier();
00507   }
00508   
00509   for (int iproc=0; iproc < NumProc; iproc++) {
00510     if (MyPID==iproc) {
00511       if (MyPID==0) {
00512   os.width(8);
00513   os <<  "   Processor ";
00514   os.width(10);
00515   os <<  "   Row Index ";
00516   os.width(10);
00517   os <<  "   Col Index ";
00518   os.width(20);
00519   os <<  "   Value     ";
00520   os << endl;
00521       }
00522       Epetra_SerialDenseVector Values(MaxNumEntries());
00523       Epetra_IntSerialDenseVector Indices(MaxNumEntries());
00524       int NumEntries;
00525       
00526       for(int i = 0; i < NumMyRows_; i++) {
00527   ExtractMyRowCopy(i, MaxNumEntries(), NumEntries, Values.Values(), Indices.Values());
00528   int Row = RowMatrixRowMap().GID(i);; // Get global row number
00529   
00530   for (int j = 0; j < NumEntries ; j++) {   
00531     int Index = RowMatrixColMap().GID(Indices[j]);
00532     os.width(8);
00533     os <<  MyPID ; os << "    ";  
00534     os.width(10);
00535     os <<  Row ; os << "    ";  
00536     os.width(10);
00537     os <<  Index; os << "    ";
00538     os.width(20);
00539     os <<  Values[j]; os << "    ";
00540     os << endl;
00541   }
00542       }
00543     
00544       os << flush;
00545       
00546     }
00547     // Do a few global ops to give I/O a chance to complete
00548     Comm().Barrier();
00549     Comm().Barrier();
00550     Comm().Barrier();
00551   }
00552   
00553   return;
00554 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines