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::ReplaceGlobalValues(int numIndices, const int* indices,
00424               const double* const* values,
00425               int format)
00426 {
00427   return(InputGlobalValues(numIndices, indices,
00428                            numIndices, indices,
00429                            values, format, REPLACE));
00430 }
00431 
00432 //----------------------------------------------------------------------------
00433 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numRows, const int* rows,
00434               int numCols, const int* cols,
00435               const double* const* values,
00436               int format)
00437 {
00438   return(InputGlobalValues(numRows, rows,
00439                            numCols, cols,
00440                            values, format, REPLACE));
00441 }
00442 
00443 //----------------------------------------------------------------------------
00444 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numIndices, const int* indices,
00445               const double* values,
00446               int format)
00447 {
00448   return(InputGlobalValues(numIndices, indices,
00449                            numIndices, indices,
00450                            values, format, REPLACE));
00451 }
00452 
00453 //----------------------------------------------------------------------------
00454 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numRows, const int* rows,
00455               int numCols, const int* cols,
00456               const double* values,
00457               int format)
00458 {
00459   return(InputGlobalValues(numRows, rows,
00460                            numCols, cols,
00461                            values, format, REPLACE));
00462 }
00463 
00464 //----------------------------------------------------------------------------
00465 int Epetra_FECrsMatrix::GlobalAssemble(bool callFillComplete)
00466 {
00467   if (Map().Comm().NumProc() < 2 || ignoreNonLocalEntries_) {
00468     if (callFillComplete) {
00469       EPETRA_CHK_ERR( FillComplete(DomainMap(), RangeMap()) );
00470     }
00471     return(0);
00472   }
00473 
00474   int i, j;
00475 
00476   //In this method we need to gather all the non-local (overlapping) data
00477   //that's been input on each processor, into the
00478   //non-overlapping distribution defined by the map that 'this' matrix was
00479   //constructed with.
00480 
00481   //First build a map that describes our nonlocal data.
00482   //We'll use the arbitrary distribution constructor of Map.
00483 
00484   Epetra_Map* sourceMap = new Epetra_Map(-1, numNonlocalRows_, nonlocalRows_,
00485            Map().IndexBase(), Map().Comm());
00486 
00487   //If sourceMap has global size 0, then no nonlocal data exists and we can
00488   //skip most of this function.
00489   if (sourceMap->NumGlobalElements() < 1) {
00490     if (callFillComplete) {
00491       EPETRA_CHK_ERR( FillComplete(DomainMap(), RangeMap()) );
00492     }
00493     delete sourceMap;
00494     return(0);
00495   }
00496 
00497   //We also need to build a column-map, containing the columns in our
00498   //nonlocal data. To do that, create a list of all column-indices that
00499   //occur in our nonlocal rows.
00500 
00501   int numCols = 0, allocLen = 0;
00502   int* cols = NULL;
00503   int insertPoint = -1;
00504 
00505   for(i=0; i<numNonlocalRows_; ++i) {
00506     for(j=0; j<nonlocalRowLengths_[i]; ++j) {
00507       int col = nonlocalCols_[i][j];
00508       int offset = Epetra_Util_binary_search(col, cols, numCols, insertPoint);
00509       if (offset < 0) {
00510   EPETRA_CHK_ERR( Epetra_Util_insert(col, insertPoint, cols,
00511              numCols, allocLen) );
00512       }
00513     }
00514   }
00515 
00516   Epetra_Map* colMap = new Epetra_Map(-1, numCols, cols,
00517               Map().IndexBase(), Map().Comm());
00518 
00519   delete [] cols;
00520   numCols = 0;
00521   allocLen = 0;
00522 
00523   //now we need to create a matrix with sourceMap and colMap, and fill it with
00524   //our nonlocal data so we can then export it to the correct owning processors.
00525 
00526   Epetra_CrsMatrix* tempMat = new Epetra_CrsMatrix(Copy, *sourceMap, *colMap,
00527                nonlocalRowLengths_);
00528 
00529 
00530   //Next we need to make sure the 'indices-are-global' attribute of tempMat's
00531   //graph is set to true, in case this processor doesn't end up calling the
00532   //InsertGlobalValues method...
00533 
00534   const Epetra_CrsGraph& graph = tempMat->Graph();
00535   Epetra_CrsGraph& nonconst_graph = const_cast<Epetra_CrsGraph&>(graph);
00536   nonconst_graph.SetIndicesAreGlobal(true);
00537 
00538   for(i=0; i<numNonlocalRows_; ++i) {
00539     EPETRA_CHK_ERR( tempMat->InsertGlobalValues(nonlocalRows_[i],
00540                  nonlocalRowLengths_[i],
00541                  nonlocalCoefs_[i],
00542                  nonlocalCols_[i]) );
00543   }
00544 
00545   //Now we need to call FillComplete on our temp matrix. We need to
00546   //pass a DomainMap and RangeMap, which are not the same as the RowMap
00547   //and ColMap that we constructed the matrix with.
00548   EPETRA_CHK_ERR(tempMat->FillComplete(DomainMap(), RangeMap()));
00549 
00550   Epetra_Export* exporter = new Epetra_Export(*sourceMap, RowMap());
00551 
00552   EPETRA_CHK_ERR(Export(*tempMat, *exporter, Add));
00553 
00554   if(callFillComplete) {
00555     EPETRA_CHK_ERR(FillComplete(DomainMap(), RangeMap()));
00556   }
00557 
00558   //now reset the values in our nonlocal data
00559   for(i=0; i<numNonlocalRows_; ++i) {
00560     for(j=0; j<nonlocalRowLengths_[i]; ++j) {
00561       nonlocalCols_[i][j] = 0;
00562       nonlocalCoefs_[i][j] = 0.0;
00563     }
00564     nonlocalRowLengths_[i] = 0;
00565   }
00566 
00567   delete exporter;
00568   delete tempMat;
00569   delete sourceMap;
00570   delete colMap;
00571 
00572   return(0);
00573 }
00574 
00575 //----------------------------------------------------------------------------
00576 int Epetra_FECrsMatrix::GlobalAssemble(const Epetra_Map& domain_map,
00577                                        const Epetra_Map& range_map,
00578                                        bool callFillComplete)
00579 {
00580   if (Map().Comm().NumProc() < 2 || ignoreNonLocalEntries_) {
00581     if (callFillComplete) {
00582       EPETRA_CHK_ERR( FillComplete(domain_map, range_map) );
00583     }
00584     return(0);
00585   }
00586 
00587   int i, j;
00588 
00589   //In this method we need to gather all the non-local (overlapping) data
00590   //that's been input on each processor, into the
00591   //non-overlapping distribution defined by the map that 'this' matrix was
00592   //constructed with.
00593 
00594   //First build a map that describes our nonlocal data.
00595   //We'll use the arbitrary distribution constructor of Map.
00596 
00597   Epetra_Map* sourceMap = new Epetra_Map(-1, numNonlocalRows_, nonlocalRows_,
00598            Map().IndexBase(), Map().Comm());
00599 
00600   //If sourceMap has global size 0, then no nonlocal data exists and we can
00601   //skip most of this function.
00602   if (sourceMap->NumGlobalElements() < 1) {
00603     if (callFillComplete) {
00604       EPETRA_CHK_ERR( FillComplete(domain_map, range_map) );
00605     }
00606     delete sourceMap;
00607     return(0);
00608   }
00609 
00610   //We also need to build a column-map, containing the columns in our
00611   //nonlocal data. To do that, create a list of all column-indices that
00612   //occur in our nonlocal rows.
00613 
00614   int numCols = 0, allocLen = 0;
00615   int* cols = NULL;
00616   int insertPoint = -1;
00617 
00618   for(i=0; i<numNonlocalRows_; ++i) {
00619     for(j=0; j<nonlocalRowLengths_[i]; ++j) {
00620       int col = nonlocalCols_[i][j];
00621       int offset = Epetra_Util_binary_search(col, cols, numCols, insertPoint);
00622       if (offset < 0) {
00623   EPETRA_CHK_ERR( Epetra_Util_insert(col, insertPoint, cols,
00624              numCols, allocLen) );
00625       }
00626     }
00627   }
00628 
00629   Epetra_Map* colMap = new Epetra_Map(-1, numCols, cols,
00630               Map().IndexBase(), Map().Comm());
00631 
00632   delete [] cols;
00633   numCols = 0;
00634   allocLen = 0;
00635 
00636   //now we need to create a matrix with sourceMap and colMap, and fill it with
00637   //our nonlocal data so we can then export it to the correct owning processors.
00638 
00639   Epetra_CrsMatrix* tempMat = new Epetra_CrsMatrix(Copy, *sourceMap, *colMap,
00640                nonlocalRowLengths_);
00641 
00642 
00643   //Next we need to make sure the 'indices-are-global' attribute of tempMat's
00644   //graph is set to true, in case this processor doesn't end up calling the
00645   //InsertGlobalValues method...
00646 
00647   const Epetra_CrsGraph& graph = tempMat->Graph();
00648   Epetra_CrsGraph& nonconst_graph = const_cast<Epetra_CrsGraph&>(graph);
00649   nonconst_graph.SetIndicesAreGlobal(true);
00650 
00651   for(i=0; i<numNonlocalRows_; ++i) {
00652     EPETRA_CHK_ERR( tempMat->InsertGlobalValues(nonlocalRows_[i],
00653                  nonlocalRowLengths_[i],
00654                  nonlocalCoefs_[i],
00655                  nonlocalCols_[i]) );
00656   }
00657 
00658   //Now we need to call FillComplete on our temp matrix. We need to
00659   //pass a DomainMap and RangeMap, which are not the same as the RowMap
00660   //and ColMap that we constructed the matrix with.
00661   EPETRA_CHK_ERR(tempMat->FillComplete(domain_map, range_map));
00662 
00663   Epetra_Export* exporter = new Epetra_Export(*sourceMap, RowMap());
00664 
00665   EPETRA_CHK_ERR(Export(*tempMat, *exporter, Add));
00666 
00667   if(callFillComplete) {
00668     EPETRA_CHK_ERR(FillComplete(domain_map, range_map));
00669   }
00670 
00671   //now reset the values in our nonlocal data
00672   for(i=0; i<numNonlocalRows_; ++i) {
00673     for(j=0; j<nonlocalRowLengths_[i]; ++j) {
00674       nonlocalCols_[i][j] = 0;
00675       nonlocalCoefs_[i][j] = 0.0;
00676     }
00677     nonlocalRowLengths_[i] = 0;
00678   }
00679 
00680   delete exporter;
00681   delete tempMat;
00682   delete sourceMap;
00683   delete colMap;
00684 
00685   return(0);
00686 }
00687 
00688 //----------------------------------------------------------------------------
00689 int Epetra_FECrsMatrix::InputGlobalValues(int numRows, const int* rows,
00690             int numCols, const int* cols,
00691             const double*const* values,
00692             int format, int mode)
00693 {
00694   if (format != Epetra_FECrsMatrix::ROW_MAJOR &&
00695       format != Epetra_FECrsMatrix::COLUMN_MAJOR) {
00696     cerr << "Epetra_FECrsMatrix: unrecognized format specifier."<<endl;
00697     return(-1);
00698   }
00699 
00700   Epetra_CrsMatrix* crsthis = static_cast<Epetra_CrsMatrix*>(this);
00701   if (crsthis == NULL) {
00702     cerr << "Epetra_FECrsMatrix::InputGlobalValues: static_cast failed."<<endl;
00703     return(-1);
00704   }
00705 
00706   if (format == Epetra_FECrsMatrix::COLUMN_MAJOR) {
00707     if (numCols > workDataLength_) {
00708       delete [] workData_;
00709       workDataLength_ = numCols*2;
00710       workData_ = new double[workDataLength_];
00711     }
00712   }
00713 
00714   int returncode = 0;
00715   int err = 0;
00716 
00717   for(int i=0; i<numRows; ++i) {
00718     double* valuesptr = (double*)values[i];
00719 
00720     if (format == Epetra_FECrsMatrix::COLUMN_MAJOR) {
00721       //if the data is in column-major order, then we need to copy the i-th
00722       //column of the values table into workData_, and that will be the i-th
00723       //row. ... Is that clear?
00724 
00725       for(int j=0; j<numCols; ++j) {
00726   workData_[j] = values[j][i];
00727       }
00728       valuesptr = workData_;
00729     }
00730 
00731     if (Map().MyGID(rows[i])) {
00732       switch(mode) {
00733       case Epetra_FECrsMatrix::SUMINTO:
00734         err = crsthis->SumIntoGlobalValues(rows[i], numCols,
00735              valuesptr, (int*)cols);
00736   if (err<0) return(err);
00737   if (err>0) returncode = err;
00738   break;
00739       case Epetra_FECrsMatrix::REPLACE:
00740   err = crsthis->ReplaceGlobalValues(rows[i], numCols,
00741              valuesptr, (int*)cols);
00742   if (err<0) return(err);
00743   if (err>0) returncode = err;
00744   break;
00745       case Epetra_FECrsMatrix::INSERT:
00746   err = crsthis->InsertGlobalValues(rows[i], numCols,
00747             valuesptr, (int*)cols);
00748   if (err<0) return(err);
00749   if (err>0) returncode = err;
00750   break;
00751       default:
00752   cerr << "Epetra_FECrsMatrix: internal error, bad input mode."<<endl;
00753   return(-1);
00754       }
00755     }
00756     else {
00757       err = InputNonlocalGlobalValues(rows[i],
00758               numCols, cols,
00759               valuesptr, mode);
00760       if (err<0) return(err);
00761       if (err>0) returncode = err;
00762     }
00763   }
00764 
00765   return(returncode);
00766 }
00767 
00768 //----------------------------------------------------------------------------
00769 int Epetra_FECrsMatrix::InputGlobalValues(int numRows, const int* rows,
00770             int numCols, const int* cols,
00771             const double* values,
00772             int format, int mode)
00773 {
00774   int first_dim = format==COLUMN_MAJOR ? numCols : numRows;
00775   int second_dim = format==COLUMN_MAJOR ? numRows : numCols;
00776 
00777   const double** values_2d = new const double*[first_dim];
00778   int offset = 0;
00779   for(int i=0; i<first_dim; ++i) {
00780     values_2d[i] = &(values[offset]);
00781     offset += second_dim;
00782   }
00783 
00784   int err = InputGlobalValues(numRows, rows, numCols, cols,
00785             values_2d, format, mode);
00786   delete [] values_2d;
00787 
00788   return(err);
00789 }
00790 
00791 //----------------------------------------------------------------------------
00792 int Epetra_FECrsMatrix::InputNonlocalGlobalValues(int row,
00793               int numCols, const int* cols,
00794               const double* values,
00795               int mode)
00796 {
00797   int insertPoint = -1;
00798 
00799   //find offset of this row in our list of nonlocal rows.
00800   int rowoffset = Epetra_Util_binary_search(row, nonlocalRows_, numNonlocalRows_,
00801               insertPoint);
00802 
00803   if (rowoffset < 0) {
00804     EPETRA_CHK_ERR( InsertNonlocalRow(row, insertPoint) );
00805     rowoffset = insertPoint;
00806   }
00807 
00808   for(int i=0; i<numCols; ++i) {
00809     EPETRA_CHK_ERR( InputNonlocalValue(rowoffset, cols[i], values[i],
00810                mode) );
00811   }
00812 
00813   return(0);
00814 }
00815 
00816 //----------------------------------------------------------------------------
00817 int Epetra_FECrsMatrix::InsertNonlocalRow(int row, int offset)
00818 {
00819   int alloc_len = numNonlocalRows_;
00820   EPETRA_CHK_ERR( Epetra_Util_insert(row, offset, nonlocalRows_, numNonlocalRows_,
00821              alloc_len, 1) );
00822 
00823   int tmp1 = numNonlocalRows_-1;
00824   int tmp2 = alloc_len-1;
00825 
00826   EPETRA_CHK_ERR( Epetra_Util_insert(0, offset, nonlocalRowLengths_,
00827              tmp1, tmp2, 1) );
00828 
00829   --tmp1;
00830   --tmp2;
00831   int initialAllocLen = 16;
00832   EPETRA_CHK_ERR( Epetra_Util_insert(initialAllocLen, offset,
00833              nonlocalRowAllocLengths_,
00834              tmp1, tmp2, 1) );
00835 
00836   int** newCols = new int*[numNonlocalRows_];
00837   double** newCoefs = new double*[numNonlocalRows_];
00838 
00839   if (newCols == NULL || newCoefs == NULL) {
00840     return(-1);
00841   }
00842 
00843   newCols[offset] = new int[initialAllocLen];
00844   newCoefs[offset] = new double[initialAllocLen];
00845 
00846   int index = 0;
00847   for(int i=0; i<numNonlocalRows_-1; ++i) {
00848     if (i == offset) {
00849       ++index;
00850     }
00851 
00852     newCols[index] = nonlocalCols_[i];
00853     newCoefs[index++] = nonlocalCoefs_[i];
00854   }
00855 
00856   delete [] nonlocalCols_;
00857   delete [] nonlocalCoefs_;
00858 
00859   nonlocalCols_ = newCols;
00860   nonlocalCoefs_ = newCoefs;
00861 
00862   return(0);
00863 }
00864 
00865 //----------------------------------------------------------------------------
00866 int Epetra_FECrsMatrix::InputNonlocalValue(int rowoffset,
00867              int col, double value,
00868              int mode)
00869 {
00870   int*& colIndices = nonlocalCols_[rowoffset];
00871   double*& coefs = nonlocalCoefs_[rowoffset];
00872 
00873   int insertPoint = -1;
00874   int coloffset = Epetra_Util_binary_search(col, colIndices,
00875               nonlocalRowLengths_[rowoffset],
00876               insertPoint);
00877 
00878   if (coloffset >= 0) {
00879     if (mode == SUMINTO || mode == INSERT) {
00880       coefs[coloffset] += value;
00881     }
00882     else {
00883       coefs[coloffset] = value;
00884     }
00885   }
00886   else {
00887     //else
00888     //  insert col in colIndices
00889     //  insert value in coefs
00890 
00891     int tmp1 = nonlocalRowLengths_[rowoffset];
00892     int tmp2 = nonlocalRowAllocLengths_[rowoffset];
00893     EPETRA_CHK_ERR( Epetra_Util_insert(col, insertPoint, colIndices, tmp1, tmp2));
00894     EPETRA_CHK_ERR( Epetra_Util_insert(value, insertPoint, coefs,
00895                nonlocalRowLengths_[rowoffset],
00896                nonlocalRowAllocLengths_[rowoffset]));
00897   }
00898 
00899   return(0);
00900 }

Generated on Thu Sep 18 12:37:57 2008 for Epetra Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1