Epetra_FECrsMatrix.cpp

Go to the documentation of this file.
00001 
00002 //@HEADER
00003 // ************************************************************************
00004 // 
00005 //               Epetra: Linear Algebra Services Package 
00006 //                 Copyright (2001) Sandia Corporation
00007 // 
00008 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00009 // license for use of this work by or on behalf of the U.S. Government.
00010 // 
00011 // This library is free software; you can redistribute it and/or modify
00012 // it under the terms of the GNU Lesser General Public License as
00013 // published by the Free Software Foundation; either version 2.1 of the
00014 // License, or (at your option) any later version.
00015 //  
00016 // This library is distributed in the hope that it will be useful, but
00017 // WITHOUT ANY WARRANTY; without even the implied warranty of
00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00019 // Lesser General Public License for more details.
00020 //  
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License along with this library; if not, write to the Free Software
00023 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00024 // USA
00025 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00026 // 
00027 // ************************************************************************
00028 //@HEADER
00029 
00030 #include <Epetra_FECrsMatrix.h>
00031 #include <Epetra_IntSerialDenseVector.h>
00032 #include <Epetra_SerialDenseMatrix.h>
00033 #include <Epetra_Export.h>
00034 #include <Epetra_Comm.h>
00035 #include <Epetra_Map.h>
00036 #include <Epetra_Util.h>
00037 
00038 //----------------------------------------------------------------------------
00039 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00040                const Epetra_Map& RowMap,
00041                int* NumEntriesPerRow,
00042                bool ignoreNonLocalEntries)
00043   : Epetra_CrsMatrix(CV, RowMap, NumEntriesPerRow),
00044     myFirstRow_(0),
00045     myNumRows_(0),
00046     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00047     numNonlocalRows_(0),
00048     nonlocalRows_(NULL),
00049     nonlocalRowLengths_(NULL),
00050     nonlocalRowAllocLengths_(NULL),
00051     nonlocalCols_(NULL),
00052     nonlocalCoefs_(NULL),
00053     workData_(NULL),
00054     workDataLength_(0)
00055 {
00056   myFirstRow_ = RowMap.MinMyGID();
00057   myNumRows_ = RowMap.NumMyElements();
00058 
00059   workData_ = new double[128];
00060   workDataLength_ = 128;
00061 }
00062 
00063 //----------------------------------------------------------------------------
00064 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00065                const Epetra_Map& RowMap,
00066                int NumEntriesPerRow,
00067                bool ignoreNonLocalEntries)
00068   : Epetra_CrsMatrix(CV, RowMap, NumEntriesPerRow),
00069     myFirstRow_(0),
00070     myNumRows_(0),
00071     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00072     numNonlocalRows_(0),
00073     nonlocalRows_(NULL),
00074     nonlocalRowLengths_(NULL),
00075     nonlocalRowAllocLengths_(NULL),
00076     nonlocalCols_(NULL),
00077     nonlocalCoefs_(NULL),
00078     workData_(NULL),
00079     workDataLength_(0)
00080 {
00081   myFirstRow_ = RowMap.MinMyGID();
00082   myNumRows_ = RowMap.NumMyElements();
00083 
00084   workData_ = new double[128];
00085   workDataLength_ = 128;
00086 }
00087 
00088 //----------------------------------------------------------------------------
00089 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00090                const Epetra_Map& RowMap,
00091                const Epetra_Map& ColMap,
00092                int* NumEntriesPerRow,
00093                bool ignoreNonLocalEntries)
00094   : Epetra_CrsMatrix(CV, RowMap, ColMap, NumEntriesPerRow),
00095     myFirstRow_(0),
00096     myNumRows_(0),
00097     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00098     numNonlocalRows_(0),
00099     nonlocalRows_(NULL),
00100     nonlocalRowLengths_(NULL),
00101     nonlocalRowAllocLengths_(NULL),
00102     nonlocalCols_(NULL),
00103     nonlocalCoefs_(NULL),
00104     workData_(NULL),
00105     workDataLength_(0)
00106 {
00107   myFirstRow_ = RowMap.MinMyGID();
00108   myNumRows_ = RowMap.NumMyElements();
00109 
00110   workData_ = new double[128];
00111   workDataLength_ = 128;
00112 }
00113 
00114 //----------------------------------------------------------------------------
00115 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00116                const Epetra_Map& RowMap,
00117                const Epetra_Map& ColMap,
00118                int NumEntriesPerRow,
00119                bool ignoreNonLocalEntries)
00120   : Epetra_CrsMatrix(CV, RowMap, ColMap, NumEntriesPerRow),
00121     myFirstRow_(0),
00122     myNumRows_(0),
00123     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00124     numNonlocalRows_(0),
00125     nonlocalRows_(NULL),
00126     nonlocalRowLengths_(NULL),
00127     nonlocalRowAllocLengths_(NULL),
00128     nonlocalCols_(NULL),
00129     nonlocalCoefs_(NULL),
00130     workData_(NULL),
00131     workDataLength_(0)
00132 {
00133   myFirstRow_ = RowMap.MinMyGID();
00134   myNumRows_ = RowMap.NumMyElements();
00135 
00136   workData_ = new double[128];
00137   workDataLength_ = 128;
00138 }
00139 
00140 //----------------------------------------------------------------------------
00141 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00142                const Epetra_CrsGraph& Graph,
00143                bool ignoreNonLocalEntries)
00144   : Epetra_CrsMatrix(CV, Graph),
00145     myFirstRow_(0),
00146     myNumRows_(0),
00147     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00148     numNonlocalRows_(0),
00149     nonlocalRows_(NULL),
00150     nonlocalRowLengths_(NULL),
00151     nonlocalRowAllocLengths_(NULL),
00152     nonlocalCols_(NULL),
00153     nonlocalCoefs_(NULL),
00154     workData_(NULL),
00155     workDataLength_(0)
00156 {
00157   myFirstRow_ = RowMap().MinMyGID();
00158   myNumRows_ = RowMap().NumMyElements();
00159 
00160   workData_ = new double[128];
00161   workDataLength_ = 128;
00162 }
00163    
00164 //----------------------------------------------------------------------------
00165 Epetra_FECrsMatrix::Epetra_FECrsMatrix(const Epetra_FECrsMatrix& src)
00166  : Epetra_CrsMatrix(src),
00167    myFirstRow_(0),
00168    myNumRows_(0),
00169    ignoreNonLocalEntries_(false),
00170    numNonlocalRows_(0),
00171    nonlocalRows_(NULL),
00172    nonlocalRowLengths_(NULL),
00173    nonlocalRowAllocLengths_(NULL),
00174    nonlocalCols_(NULL),
00175    nonlocalCoefs_(NULL),
00176    workData_(NULL),
00177    workDataLength_(0)
00178 {
00179   operator=(src);
00180 }
00181 
00182 //----------------------------------------------------------------------------
00183 Epetra_FECrsMatrix& Epetra_FECrsMatrix::operator=(const Epetra_FECrsMatrix& src)
00184 {
00185   if (this == &src) {
00186     return( *this );
00187   }
00188 
00189   DeleteMemory();
00190 
00191   Epetra_CrsMatrix::operator=(src);
00192 
00193   myFirstRow_ = src.myFirstRow_;
00194   myNumRows_ = src.myNumRows_;
00195   ignoreNonLocalEntries_ = src.ignoreNonLocalEntries_;
00196   numNonlocalRows_ = src.numNonlocalRows_;
00197 
00198   workDataLength_ = 128;
00199   workData_ = new double[workDataLength_];
00200 
00201   if (numNonlocalRows_ < 1) {
00202     return( *this );
00203   }
00204 
00205   nonlocalRows_ = new int[numNonlocalRows_];
00206   nonlocalRowLengths_ = new int[numNonlocalRows_];
00207   nonlocalRowAllocLengths_ = new int[numNonlocalRows_];
00208   nonlocalCols_ = new int*[numNonlocalRows_];
00209   nonlocalCoefs_= new double*[numNonlocalRows_];
00210 
00211   for(int i=0; i<numNonlocalRows_; ++i) {
00212     nonlocalRows_[i] = src.nonlocalRows_[i];
00213     nonlocalRowLengths_[i] = src.nonlocalRowLengths_[i];
00214     nonlocalRowAllocLengths_[i] = src.nonlocalRowAllocLengths_[i];
00215 
00216     nonlocalCols_[i] = new int[nonlocalRowAllocLengths_[i]];
00217     nonlocalCoefs_[i] = new double[nonlocalRowAllocLengths_[i]];
00218 
00219     for(int j=0; j<nonlocalRowLengths_[i]; ++j) {
00220       nonlocalCols_[i][j] = src.nonlocalCols_[i][j];
00221       nonlocalCoefs_[i][j] = src.nonlocalCoefs_[i][j];
00222     }
00223   }
00224 
00225   return( *this );
00226 }
00227 
00228 //----------------------------------------------------------------------------
00229 Epetra_FECrsMatrix::~Epetra_FECrsMatrix()
00230 {
00231   DeleteMemory();
00232 }
00233 
00234 //----------------------------------------------------------------------------
00235 void Epetra_FECrsMatrix::DeleteMemory()
00236 {
00237   if (numNonlocalRows_ > 0) {
00238     for(int i=0; i<numNonlocalRows_; ++i) {
00239       delete [] nonlocalCols_[i];
00240       delete [] nonlocalCoefs_[i];
00241     }
00242     delete [] nonlocalCols_;
00243     delete [] nonlocalCoefs_;
00244     delete [] nonlocalRows_;
00245     delete [] nonlocalRowLengths_;
00246     delete [] nonlocalRowAllocLengths_;
00247     numNonlocalRows_ = 0;
00248   }
00249 
00250   delete [] workData_;
00251   workDataLength_ = 0;
00252 }
00253 
00254 //----------------------------------------------------------------------------
00255 int Epetra_FECrsMatrix::SumIntoGlobalValues(int numIndices, const int* indices,
00256               const double* const* values,
00257               int format)
00258 {
00259   return(InputGlobalValues(numIndices, indices,
00260                            numIndices, indices,
00261                            values, format, SUMINTO));
00262 }
00263 
00264 //----------------------------------------------------------------------------
00265 int Epetra_FECrsMatrix::SumIntoGlobalValues(int numRows, const int* rows,
00266               int numCols, const int* cols,
00267               const double* const* values,
00268               int format)
00269 {
00270   return(InputGlobalValues(numRows, rows,
00271                            numCols, cols,
00272                            values, format, SUMINTO));
00273 }
00274 
00275 //----------------------------------------------------------------------------
00276 int Epetra_FECrsMatrix::SumIntoGlobalValues(int numIndices, const int* indices,
00277               const double* values,
00278               int format)
00279 {
00280   return(InputGlobalValues(numIndices, indices,
00281                            numIndices, indices,
00282                            values, format, SUMINTO));
00283 }
00284 
00285 //----------------------------------------------------------------------------
00286 int Epetra_FECrsMatrix::SumIntoGlobalValues(int numRows, const int* rows,
00287               int numCols, const int* cols,
00288               const double* values,
00289               int format)
00290 {
00291   return(InputGlobalValues(numRows, rows,
00292                            numCols, cols,
00293                            values, format, SUMINTO));
00294 }
00295 
00296 //----------------------------------------------------------------------------
00297 int Epetra_FECrsMatrix::SumIntoGlobalValues(const Epetra_IntSerialDenseVector& indices,
00298               const Epetra_SerialDenseMatrix& values,
00299               int format)
00300 {
00301   if (indices.Length() != values.M() || indices.Length() != values.N()) {
00302     return(-1);
00303   }
00304 
00305   return( SumIntoGlobalValues(indices.Length(), indices.Values(),
00306             values.A(), format) );
00307 }
00308 
00309 //----------------------------------------------------------------------------
00310 int Epetra_FECrsMatrix::InsertGlobalValues(const Epetra_IntSerialDenseVector& indices,
00311               const Epetra_SerialDenseMatrix& values,
00312               int format)
00313 {
00314   if (indices.Length() != values.M() || indices.Length() != values.N()) {
00315     return(-1);
00316   }
00317 
00318   return( InsertGlobalValues(indices.Length(), indices.Values(),
00319             values.A(), format) );
00320 }
00321 
00322 //----------------------------------------------------------------------------
00323 int Epetra_FECrsMatrix::ReplaceGlobalValues(const Epetra_IntSerialDenseVector& indices,
00324               const Epetra_SerialDenseMatrix& values,
00325               int format)
00326 {
00327   if (indices.Length() != values.M() || indices.Length() != values.N()) {
00328     return(-1);
00329   }
00330 
00331   return( ReplaceGlobalValues(indices.Length(), indices.Values(),
00332             values.A(), format) );
00333 }
00334 
00335 //----------------------------------------------------------------------------
00336 int Epetra_FECrsMatrix::SumIntoGlobalValues(const Epetra_IntSerialDenseVector& rows,
00337               const Epetra_IntSerialDenseVector& cols,
00338               const Epetra_SerialDenseMatrix& values,
00339               int format)
00340 {
00341   if (rows.Length() != values.M() || cols.Length() != values.N()) {
00342     return(-1);
00343   }
00344 
00345   return( SumIntoGlobalValues(rows.Length(), rows.Values(),
00346             cols.Length(), cols.Values(),
00347             values.A(), format) );
00348 }
00349 
00350 //----------------------------------------------------------------------------
00351 int Epetra_FECrsMatrix::InsertGlobalValues(const Epetra_IntSerialDenseVector& rows,
00352              const Epetra_IntSerialDenseVector& cols,
00353              const Epetra_SerialDenseMatrix& values,
00354              int format)
00355 {
00356   if (rows.Length() != values.M() || cols.Length() != values.N()) {
00357     return(-1);
00358   }
00359 
00360   return( InsertGlobalValues(rows.Length(), rows.Values(),
00361            cols.Length(), cols.Values(),
00362             values.A(), format) );
00363 }
00364 
00365 //----------------------------------------------------------------------------
00366 int Epetra_FECrsMatrix::ReplaceGlobalValues(const Epetra_IntSerialDenseVector& rows,
00367               const Epetra_IntSerialDenseVector& cols,
00368               const Epetra_SerialDenseMatrix& values,
00369               int format)
00370 {
00371   if (rows.Length() != values.M() || cols.Length() != values.N()) {
00372     return(-1);
00373   }
00374 
00375   return( ReplaceGlobalValues(rows.Length(), rows.Values(),
00376             cols.Length(), cols.Values(),
00377             values.A(), format) );
00378 }
00379 
00380 //----------------------------------------------------------------------------
00381 int Epetra_FECrsMatrix::InsertGlobalValues(int numIndices, const int* indices,
00382               const double* const* values,
00383               int format)
00384 {
00385   return(InputGlobalValues(numIndices, indices,
00386                            numIndices, indices,
00387                            values, format, INSERT));
00388 }
00389 
00390 //----------------------------------------------------------------------------
00391 int Epetra_FECrsMatrix::InsertGlobalValues(int numRows, const int* rows,
00392               int numCols, const int* cols,
00393               const double* const* values,
00394               int format)
00395 {
00396   return(InputGlobalValues(numRows, rows,
00397                            numCols, cols,
00398                            values, format, INSERT));
00399 }
00400 
00401 //----------------------------------------------------------------------------
00402 int Epetra_FECrsMatrix::InsertGlobalValues(int numIndices, const int* indices,
00403               const double* values,
00404               int format)
00405 {
00406   return(InputGlobalValues(numIndices, indices,
00407                            numIndices, indices,
00408                            values, format, INSERT));
00409 }
00410 
00411 //----------------------------------------------------------------------------
00412 int Epetra_FECrsMatrix::InsertGlobalValues(int numRows, const int* rows,
00413               int numCols, const int* cols,
00414               const double* values,
00415               int format)
00416 {
00417   return(InputGlobalValues(numRows, rows,
00418                            numCols, cols,
00419                            values, format, INSERT));
00420 }
00421 
00422 //----------------------------------------------------------------------------
00423 int Epetra_FECrsMatrix::SumIntoGlobalValues(int GlobalRow, int NumEntries,
00424                                             double* Values, int* Indices)
00425 {
00426   return(InputGlobalValues(1, &GlobalRow,
00427                            NumEntries, Indices, Values,
00428                            Epetra_FECrsMatrix::ROW_MAJOR, SUMINTO));
00429 }
00430 
00431 //----------------------------------------------------------------------------
00432 int Epetra_FECrsMatrix::InsertGlobalValues(int GlobalRow, int NumEntries,
00433                                             double* Values, int* Indices)
00434 {
00435   return(InputGlobalValues(1, &GlobalRow,
00436                            NumEntries, Indices, Values,
00437                            Epetra_FECrsMatrix::ROW_MAJOR, INSERT));
00438 }
00439 
00440 //----------------------------------------------------------------------------
00441 int Epetra_FECrsMatrix::ReplaceGlobalValues(int GlobalRow, int NumEntries,
00442                                             double* Values, int* Indices)
00443 {
00444   return(InputGlobalValues(1, &GlobalRow,
00445                            NumEntries, Indices, Values,
00446                            Epetra_FECrsMatrix::ROW_MAJOR, REPLACE));
00447 }
00448 
00449 //----------------------------------------------------------------------------
00450 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numIndices, const int* indices,
00451               const double* const* values,
00452               int format)
00453 {
00454   return(InputGlobalValues(numIndices, indices,
00455                            numIndices, indices,
00456                            values, format, REPLACE));
00457 }
00458 
00459 //----------------------------------------------------------------------------
00460 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numRows, const int* rows,
00461               int numCols, const int* cols,
00462               const double* const* values,
00463               int format)
00464 {
00465   return(InputGlobalValues(numRows, rows,
00466                            numCols, cols,
00467                            values, format, REPLACE));
00468 }
00469 
00470 //----------------------------------------------------------------------------
00471 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numIndices, const int* indices,
00472               const double* values,
00473               int format)
00474 {
00475   return(InputGlobalValues(numIndices, indices,
00476                            numIndices, indices,
00477                            values, format, REPLACE));
00478 }
00479 
00480 //----------------------------------------------------------------------------
00481 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numRows, const int* rows,
00482               int numCols, const int* cols,
00483               const double* values,
00484               int format)
00485 {
00486   return(InputGlobalValues(numRows, rows,
00487                            numCols, cols,
00488                            values, format, REPLACE));
00489 }
00490 
00491 //----------------------------------------------------------------------------
00492 int Epetra_FECrsMatrix::GlobalAssemble(bool callFillComplete)
00493 {
00494   if (Map().Comm().NumProc() < 2 || ignoreNonLocalEntries_) {
00495     if (callFillComplete) {
00496       EPETRA_CHK_ERR( FillComplete(DomainMap(), RangeMap()) );
00497     }
00498     return(0);
00499   }
00500 
00501   //In this method we need to gather all the non-local (overlapping) data
00502   //that's been input on each processor, into the
00503   //non-overlapping distribution defined by the map that 'this' matrix was
00504   //constructed with.
00505 
00506   //First build a map that describes our nonlocal data.
00507   //We'll use the arbitrary distribution constructor of Map.
00508 
00509   Epetra_Map sourceMap(-1, numNonlocalRows_, nonlocalRows_,
00510        Map().IndexBase(), Map().Comm());
00511 
00512   //If sourceMap has global size 0, then no nonlocal data exists and we can
00513   //skip most of this function.
00514   if (sourceMap.NumGlobalElements() < 1) {
00515     if (callFillComplete) {
00516       EPETRA_CHK_ERR( FillComplete(DomainMap(), RangeMap()) );
00517     }
00518     return(0);
00519   }
00520 
00521   //We also need to build a column-map, containing the columns in our
00522   //nonlocal data. To do that, create a list of all column-indices that
00523   //occur in our nonlocal rows.
00524 
00525   int numCols = 0, allocLen = 0;
00526   int* cols = NULL;
00527   int insertPoint = -1;
00528 
00529   for(int i=0; i<numNonlocalRows_; ++i) {
00530     for(int j=0; j<nonlocalRowLengths_[i]; ++j) {
00531       int col = nonlocalCols_[i][j];
00532       int offset = Epetra_Util_binary_search(col, cols, numCols, insertPoint);
00533       if (offset < 0) {
00534   EPETRA_CHK_ERR( Epetra_Util_insert(col, insertPoint, cols,
00535              numCols, allocLen) );
00536       }
00537     }
00538   }
00539 
00540   Epetra_Map colMap(-1, numCols, cols,
00541           Map().IndexBase(), Map().Comm());
00542 
00543   delete [] cols;
00544   numCols = 0;
00545   allocLen = 0;
00546 
00547   //now we need to create a matrix with sourceMap and colMap, and fill it with
00548   //our nonlocal data so we can then export it to the correct owning processors.
00549 
00550   Epetra_CrsMatrix tempMat(Copy, sourceMap, colMap, nonlocalRowLengths_);
00551 
00552   //Next we need to make sure the 'indices-are-global' attribute of tempMat's
00553   //graph is set to true, in case this processor doesn't end up calling the
00554   //InsertGlobalValues method...
00555 
00556   const Epetra_CrsGraph& graph = tempMat.Graph();
00557   Epetra_CrsGraph& nonconst_graph = const_cast<Epetra_CrsGraph&>(graph);
00558   nonconst_graph.SetIndicesAreGlobal(true);
00559 
00560   for(int i=0; i<numNonlocalRows_; ++i) {
00561     EPETRA_CHK_ERR( tempMat.InsertGlobalValues(nonlocalRows_[i],
00562                  nonlocalRowLengths_[i],
00563                  nonlocalCoefs_[i],
00564                  nonlocalCols_[i]) );
00565   }
00566 
00567   //Now we need to call FillComplete on our temp matrix. We need to
00568   //pass a DomainMap and RangeMap, which are not the same as the RowMap
00569   //and ColMap that we constructed the matrix with.
00570   EPETRA_CHK_ERR(tempMat.FillComplete(DomainMap(), RangeMap()));
00571 
00572   Epetra_Export exporter(sourceMap, RowMap());
00573 
00574   EPETRA_CHK_ERR(Export(tempMat, exporter, Add));
00575 
00576   if(callFillComplete) {
00577     EPETRA_CHK_ERR(FillComplete(DomainMap(), RangeMap()));
00578   }
00579 
00580   //now reset the values in our nonlocal data
00581   for(int i=0; i<numNonlocalRows_; ++i) {
00582     for(int j=0; j<nonlocalRowLengths_[i]; ++j) {
00583       nonlocalCols_[i][j] = 0;
00584       nonlocalCoefs_[i][j] = 0.0;
00585     }
00586     nonlocalRowLengths_[i] = 0;
00587   }
00588 
00589   return(0);
00590 }
00591 
00592 //----------------------------------------------------------------------------
00593 int Epetra_FECrsMatrix::GlobalAssemble(const Epetra_Map& domain_map,
00594                                        const Epetra_Map& range_map,
00595                                        bool callFillComplete)
00596 {
00597   if (Map().Comm().NumProc() < 2 || ignoreNonLocalEntries_) {
00598     if (callFillComplete) {
00599       EPETRA_CHK_ERR( FillComplete(domain_map, range_map) );
00600     }
00601     return(0);
00602   }
00603 
00604   int i, j;
00605 
00606   //In this method we need to gather all the non-local (overlapping) data
00607   //that's been input on each processor, into the
00608   //non-overlapping distribution defined by the map that 'this' matrix was
00609   //constructed with.
00610 
00611   //First build a map that describes our nonlocal data.
00612   //We'll use the arbitrary distribution constructor of Map.
00613 
00614   Epetra_Map* sourceMap = new Epetra_Map(-1, numNonlocalRows_, nonlocalRows_,
00615            Map().IndexBase(), Map().Comm());
00616 
00617   //If sourceMap has global size 0, then no nonlocal data exists and we can
00618   //skip most of this function.
00619   if (sourceMap->NumGlobalElements() < 1) {
00620     if (callFillComplete) {
00621       EPETRA_CHK_ERR( FillComplete(domain_map, range_map) );
00622     }
00623     delete sourceMap;
00624     return(0);
00625   }
00626 
00627   //We also need to build a column-map, containing the columns in our
00628   //nonlocal data. To do that, create a list of all column-indices that
00629   //occur in our nonlocal rows.
00630 
00631   int numCols = 0, allocLen = 0;
00632   int* cols = NULL;
00633   int insertPoint = -1;
00634 
00635   for(i=0; i<numNonlocalRows_; ++i) {
00636     for(j=0; j<nonlocalRowLengths_[i]; ++j) {
00637       int col = nonlocalCols_[i][j];
00638       int offset = Epetra_Util_binary_search(col, cols, numCols, insertPoint);
00639       if (offset < 0) {
00640   EPETRA_CHK_ERR( Epetra_Util_insert(col, insertPoint, cols,
00641              numCols, allocLen) );
00642       }
00643     }
00644   }
00645 
00646   Epetra_Map* colMap = new Epetra_Map(-1, numCols, cols,
00647               Map().IndexBase(), Map().Comm());
00648 
00649   delete [] cols;
00650   numCols = 0;
00651   allocLen = 0;
00652 
00653   //now we need to create a matrix with sourceMap and colMap, and fill it with
00654   //our nonlocal data so we can then export it to the correct owning processors.
00655 
00656   Epetra_CrsMatrix* tempMat = new Epetra_CrsMatrix(Copy, *sourceMap, *colMap,
00657                nonlocalRowLengths_);
00658 
00659 
00660   //Next we need to make sure the 'indices-are-global' attribute of tempMat's
00661   //graph is set to true, in case this processor doesn't end up calling the
00662   //InsertGlobalValues method...
00663 
00664   const Epetra_CrsGraph& graph = tempMat->Graph();
00665   Epetra_CrsGraph& nonconst_graph = const_cast<Epetra_CrsGraph&>(graph);
00666   nonconst_graph.SetIndicesAreGlobal(true);
00667 
00668   for(i=0; i<numNonlocalRows_; ++i) {
00669     EPETRA_CHK_ERR( tempMat->InsertGlobalValues(nonlocalRows_[i],
00670                  nonlocalRowLengths_[i],
00671                  nonlocalCoefs_[i],
00672                  nonlocalCols_[i]) );
00673   }
00674 
00675   //Now we need to call FillComplete on our temp matrix. We need to
00676   //pass a DomainMap and RangeMap, which are not the same as the RowMap
00677   //and ColMap that we constructed the matrix with.
00678   EPETRA_CHK_ERR(tempMat->FillComplete(domain_map, range_map));
00679 
00680   Epetra_Export* exporter = new Epetra_Export(*sourceMap, RowMap());
00681 
00682   EPETRA_CHK_ERR(Export(*tempMat, *exporter, Add));
00683 
00684   if(callFillComplete) {
00685     EPETRA_CHK_ERR(FillComplete(domain_map, range_map));
00686   }
00687 
00688   //now reset the values in our nonlocal data
00689   for(i=0; i<numNonlocalRows_; ++i) {
00690     for(j=0; j<nonlocalRowLengths_[i]; ++j) {
00691       nonlocalCols_[i][j] = 0;
00692       nonlocalCoefs_[i][j] = 0.0;
00693     }
00694     nonlocalRowLengths_[i] = 0;
00695   }
00696 
00697   delete exporter;
00698   delete tempMat;
00699   delete sourceMap;
00700   delete colMap;
00701 
00702   return(0);
00703 }
00704 
00705 //----------------------------------------------------------------------------
00706 int Epetra_FECrsMatrix::InputGlobalValues(int numRows, const int* rows,
00707             int numCols, const int* cols,
00708             const double*const* values,
00709             int format, int mode)
00710 {
00711   if (format != Epetra_FECrsMatrix::ROW_MAJOR &&
00712       format != Epetra_FECrsMatrix::COLUMN_MAJOR) {
00713     cerr << "Epetra_FECrsMatrix: unrecognized format specifier."<<endl;
00714     return(-1);
00715   }
00716 
00717   if (format == Epetra_FECrsMatrix::COLUMN_MAJOR) {
00718     if (numCols > workDataLength_) {
00719       delete [] workData_;
00720       workDataLength_ = numCols*2;
00721       workData_ = new double[workDataLength_];
00722     }
00723   }
00724 
00725   int returncode = 0;
00726   int err = 0;
00727 
00728   for(int i=0; i<numRows; ++i) {
00729     double* valuesptr = (double*)values[i];
00730 
00731     if (format == Epetra_FECrsMatrix::COLUMN_MAJOR) {
00732       //if the data is in column-major order, then we need to copy the i-th
00733       //column of the values table into workData_, and that will be the i-th
00734       //row. ... Is that clear?
00735 
00736       for(int j=0; j<numCols; ++j) {
00737   workData_[j] = values[j][i];
00738       }
00739       valuesptr = workData_;
00740     }
00741 
00742     if (Map().MyGID(rows[i])) {
00743       switch(mode) {
00744       case Epetra_FECrsMatrix::SUMINTO:
00745         err = this->Epetra_CrsMatrix::SumIntoGlobalValues(rows[i], numCols,
00746              valuesptr, (int*)cols);
00747   if (err<0) return(err);
00748   if (err>0) returncode = err;
00749   break;
00750       case Epetra_FECrsMatrix::REPLACE:
00751   err = this->Epetra_CrsMatrix::ReplaceGlobalValues(rows[i], numCols,
00752              valuesptr, (int*)cols);
00753   if (err<0) return(err);
00754   if (err>0) returncode = err;
00755   break;
00756       case Epetra_FECrsMatrix::INSERT:
00757   err = this->Epetra_CrsMatrix::InsertGlobalValues(rows[i], numCols,
00758             valuesptr, (int*)cols);
00759   if (err<0) return(err);
00760   if (err>0) returncode = err;
00761   break;
00762       default:
00763   cerr << "Epetra_FECrsMatrix: internal error, bad input mode."<<endl;
00764   return(-1);
00765       }
00766     }
00767     else {
00768       err = InputNonlocalGlobalValues(rows[i],
00769               numCols, cols,
00770               valuesptr, mode);
00771       if (err<0) return(err);
00772       if (err>0) returncode = err;
00773     }
00774   }
00775 
00776   return(returncode);
00777 }
00778 
00779 //----------------------------------------------------------------------------
00780 int Epetra_FECrsMatrix::InputGlobalValues(int numRows, const int* rows,
00781             int numCols, const int* cols,
00782             const double* values,
00783             int format, int mode)
00784 {
00785   int first_dim = format==COLUMN_MAJOR ? numCols : numRows;
00786   int second_dim = format==COLUMN_MAJOR ? numRows : numCols;
00787 
00788   const double** values_2d = new const double*[first_dim];
00789   int offset = 0;
00790   for(int i=0; i<first_dim; ++i) {
00791     values_2d[i] = &(values[offset]);
00792     offset += second_dim;
00793   }
00794 
00795   int err = InputGlobalValues(numRows, rows, numCols, cols,
00796             values_2d, format, mode);
00797   delete [] values_2d;
00798 
00799   return(err);
00800 }
00801 
00802 //----------------------------------------------------------------------------
00803 int Epetra_FECrsMatrix::InputNonlocalGlobalValues(int row,
00804               int numCols, const int* cols,
00805               const double* values,
00806               int mode)
00807 {
00808   int insertPoint = -1;
00809 
00810   //find offset of this row in our list of nonlocal rows.
00811   int rowoffset = Epetra_Util_binary_search(row, nonlocalRows_,
00812                                           numNonlocalRows_, insertPoint);
00813 
00814   if (rowoffset < 0) {
00815     EPETRA_CHK_ERR( InsertNonlocalRow(row, insertPoint) );
00816     rowoffset = insertPoint;
00817   }
00818 
00819   for(int i=0; i<numCols; ++i) {
00820     EPETRA_CHK_ERR( InputNonlocalValue(rowoffset, cols[i], values[i],
00821                mode) );
00822   }
00823 
00824   return(0);
00825 }
00826 
00827 //----------------------------------------------------------------------------
00828 int Epetra_FECrsMatrix::InsertNonlocalRow(int row, int offset)
00829 {
00830   int alloc_len = numNonlocalRows_;
00831   EPETRA_CHK_ERR( Epetra_Util_insert(row, offset, nonlocalRows_,
00832                                      numNonlocalRows_, alloc_len, 1) );
00833 
00834   int tmp1 = numNonlocalRows_-1;
00835   int tmp2 = alloc_len-1;
00836 
00837   EPETRA_CHK_ERR( Epetra_Util_insert(0, offset, nonlocalRowLengths_,
00838              tmp1, tmp2, 1) );
00839 
00840   --tmp1;
00841   --tmp2;
00842   int initialAllocLen = 16;
00843   EPETRA_CHK_ERR( Epetra_Util_insert(initialAllocLen, offset,
00844              nonlocalRowAllocLengths_,
00845              tmp1, tmp2, 1) );
00846 
00847   int** newCols = new int*[numNonlocalRows_];
00848   double** newCoefs = new double*[numNonlocalRows_];
00849 
00850   if (newCols == NULL || newCoefs == NULL) {
00851     return(-1);
00852   }
00853 
00854   newCols[offset] = new int[initialAllocLen];
00855   newCoefs[offset] = new double[initialAllocLen];
00856 
00857   int index = 0;
00858   for(int i=0; i<numNonlocalRows_-1; ++i) {
00859     if (i == offset) {
00860       ++index;
00861     }
00862 
00863     newCols[index] = nonlocalCols_[i];
00864     newCoefs[index++] = nonlocalCoefs_[i];
00865   }
00866 
00867   delete [] nonlocalCols_;
00868   delete [] nonlocalCoefs_;
00869 
00870   nonlocalCols_ = newCols;
00871   nonlocalCoefs_ = newCoefs;
00872 
00873   return(0);
00874 }
00875 
00876 //----------------------------------------------------------------------------
00877 int Epetra_FECrsMatrix::InputNonlocalValue(int rowoffset,
00878              int col, double value,
00879              int mode)
00880 {
00881   int*& colIndices = nonlocalCols_[rowoffset];
00882   double*& coefs = nonlocalCoefs_[rowoffset];
00883   int len = nonlocalRowLengths_[rowoffset];
00884 
00885   int insertPoint = -1;
00886   int coloffset = Epetra_Util_binary_search(col, colIndices,
00887               len, insertPoint);
00888 
00889   if (coloffset >= 0) {
00890     if (mode == SUMINTO || mode == INSERT) {
00891       coefs[coloffset] += value;
00892     }
00893     else {
00894       coefs[coloffset] = value;
00895     }
00896   }
00897   else {
00898     //else
00899     //  insert col in colIndices
00900     //  insert value in coefs
00901 
00902     int tmp1 = nonlocalRowLengths_[rowoffset];
00903     int tmp2 = nonlocalRowAllocLengths_[rowoffset];
00904     EPETRA_CHK_ERR( Epetra_Util_insert(col, insertPoint, colIndices, tmp1, tmp2));
00905     EPETRA_CHK_ERR( Epetra_Util_insert(value, insertPoint, coefs,
00906                nonlocalRowLengths_[rowoffset],
00907                nonlocalRowAllocLengths_[rowoffset]));
00908   }
00909 
00910   return(0);
00911 }
00912 

Generated on Wed May 12 21:41:05 2010 for Epetra Package Browser (Single Doxygen Collection) by  doxygen 1.4.7