Epetra Package Browser (Single Doxygen Collection) Development
Epetra_FECrsMatrix.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_FECrsMatrix.h>
00044 #include <Epetra_IntSerialDenseVector.h>
00045 #include <Epetra_SerialDenseMatrix.h>
00046 #include <Epetra_FECrsGraph.h>
00047 #include <Epetra_Export.h>
00048 #include <Epetra_Comm.h>
00049 #include <Epetra_Map.h>
00050 #include <Epetra_Util.h>
00051 
00052 //----------------------------------------------------------------------------
00053 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00054                const Epetra_Map& rowMap,
00055                int* NumEntriesPerRow,
00056                bool ignoreNonLocalEntries)
00057   : Epetra_CrsMatrix(CV, rowMap, NumEntriesPerRow),
00058     myFirstRow_(0),
00059     myNumRows_(0),
00060     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00061     nonlocalRows_(),
00062     nonlocalCols_(),
00063     nonlocalCoefs_(),
00064     workData_(128),
00065     useNonlocalMatrix_ (false),
00066     nonlocalMatrix_ (NULL)
00067 {
00068   myFirstRow_ = rowMap.MinMyGID();
00069   myNumRows_ = rowMap.NumMyElements();
00070 }
00071 
00072 //----------------------------------------------------------------------------
00073 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00074                const Epetra_Map& rowMap,
00075                int NumEntriesPerRow,
00076                bool ignoreNonLocalEntries)
00077   : Epetra_CrsMatrix(CV, rowMap, NumEntriesPerRow),
00078     myFirstRow_(0),
00079     myNumRows_(0),
00080     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00081     nonlocalRows_(),
00082     nonlocalCols_(),
00083     nonlocalCoefs_(),
00084     workData_(128),
00085     useNonlocalMatrix_ (false),
00086     nonlocalMatrix_ (NULL)
00087 {
00088   myFirstRow_ = rowMap.MinMyGID();
00089   myNumRows_ = rowMap.NumMyElements();
00090 }
00091 
00092 //----------------------------------------------------------------------------
00093 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00094                const Epetra_Map& rowMap,
00095                const Epetra_Map& colMap,
00096                int* NumEntriesPerRow,
00097                bool ignoreNonLocalEntries)
00098   : Epetra_CrsMatrix(CV, rowMap, colMap, NumEntriesPerRow),
00099     myFirstRow_(0),
00100     myNumRows_(0),
00101     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00102     nonlocalRows_(),
00103     nonlocalCols_(),
00104     nonlocalCoefs_(),
00105     workData_(128),
00106     useNonlocalMatrix_ (false),
00107     nonlocalMatrix_ (NULL)
00108 {
00109   myFirstRow_ = rowMap.MinMyGID();
00110   myNumRows_ = rowMap.NumMyElements();
00111 }
00112 
00113 //----------------------------------------------------------------------------
00114 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00115                const Epetra_Map& rowMap,
00116                const Epetra_Map& colMap,
00117                int NumEntriesPerRow,
00118                bool ignoreNonLocalEntries)
00119   : Epetra_CrsMatrix(CV, rowMap, colMap, NumEntriesPerRow),
00120     myFirstRow_(0),
00121     myNumRows_(0),
00122     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00123     nonlocalRows_(),
00124     nonlocalCols_(),
00125     nonlocalCoefs_(),
00126     workData_(128),
00127     useNonlocalMatrix_ (false),
00128     nonlocalMatrix_ (NULL)
00129 {
00130   myFirstRow_ = rowMap.MinMyGID();
00131   myNumRows_ = rowMap.NumMyElements();
00132 }
00133 
00134 //----------------------------------------------------------------------------
00135 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00136                const Epetra_CrsGraph& graph,
00137                bool ignoreNonLocalEntries)
00138   : Epetra_CrsMatrix(CV, graph),
00139     myFirstRow_(0),
00140     myNumRows_(0),
00141     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00142     nonlocalRows_(),
00143     nonlocalCols_(),
00144     nonlocalCoefs_(),
00145     workData_(128),
00146     useNonlocalMatrix_ (false),
00147     nonlocalMatrix_ (NULL)
00148 {
00149   myFirstRow_ = RowMap().MinMyGID();
00150   myNumRows_ = RowMap().NumMyElements();
00151 }
00152    
00153 //----------------------------------------------------------------------------
00154 Epetra_FECrsMatrix::Epetra_FECrsMatrix(Epetra_DataAccess CV,
00155               const Epetra_FECrsGraph& graph,
00156               bool ignoreNonLocalEntries)
00157   : Epetra_CrsMatrix(CV, graph),
00158     myFirstRow_(0),
00159     myNumRows_(0),
00160     ignoreNonLocalEntries_(ignoreNonLocalEntries),
00161     nonlocalRows_(),
00162     nonlocalCols_(),
00163     nonlocalCoefs_(),
00164     workData_(128),
00165     useNonlocalMatrix_ (graph.UseNonlocalGraph() && graph.nonlocalGraph_ != 0),
00166     nonlocalMatrix_ (useNonlocalMatrix_ ? 
00167         new Epetra_CrsMatrix(Copy,*graph.nonlocalGraph_) : NULL)
00168 {
00169   myFirstRow_ = RowMap().MinMyGID();
00170   myNumRows_ = RowMap().NumMyElements();
00171 }
00172    
00173 //----------------------------------------------------------------------------
00174 Epetra_FECrsMatrix::Epetra_FECrsMatrix(const Epetra_FECrsMatrix& src)
00175  : Epetra_CrsMatrix(src),
00176    myFirstRow_(0),
00177    myNumRows_(0),
00178    ignoreNonLocalEntries_(false),
00179    nonlocalRows_(),
00180    nonlocalCols_(),
00181    nonlocalCoefs_(),
00182    workData_(128),
00183    nonlocalMatrix_ (NULL)
00184 {
00185   operator=(src);
00186 }
00187 
00188 //----------------------------------------------------------------------------
00189 Epetra_FECrsMatrix& Epetra_FECrsMatrix::operator=(const Epetra_FECrsMatrix& src)
00190 {
00191   if (this == &src) {
00192     return( *this );
00193   }
00194 
00195   DeleteMemory();
00196 
00197   Epetra_CrsMatrix::operator=(src);
00198 
00199   useNonlocalMatrix_ = src.useNonlocalMatrix_;
00200 
00201   myFirstRow_ = src.myFirstRow_;
00202   myNumRows_ = src.myNumRows_;
00203   ignoreNonLocalEntries_ = src.ignoreNonLocalEntries_;
00204 
00205   if (nonlocalRows_.size() < 1) {
00206     return( *this );
00207   }
00208 
00209   if (useNonlocalMatrix_ && src.nonlocalMatrix_ != 0) {
00210     *nonlocalMatrix_ = *src.nonlocalMatrix_;
00211     return( *this );
00212   }
00213 
00214   nonlocalRows_ = src.nonlocalRows_;
00215   nonlocalCols_ = src.nonlocalCols_;
00216   nonlocalCoefs_= src.nonlocalCoefs_;
00217 
00218   return( *this );
00219 }
00220 
00221 //----------------------------------------------------------------------------
00222 Epetra_FECrsMatrix::~Epetra_FECrsMatrix()
00223 {
00224   DeleteMemory();
00225 }
00226 
00227 //----------------------------------------------------------------------------
00228 void Epetra_FECrsMatrix::DeleteMemory()
00229 {
00230   nonlocalRows_.clear();
00231   nonlocalCols_.clear();
00232   nonlocalCoefs_.clear();
00233 
00234   if (nonlocalMatrix_ != 0)
00235     delete nonlocalMatrix_;
00236 }
00237 
00238 //----------------------------------------------------------------------------
00239 int Epetra_FECrsMatrix::SumIntoGlobalValues(int numIndices, const int* indices,
00240               const double* const* values,
00241               int format)
00242 {
00243   return(InputGlobalValues(numIndices, indices,
00244                            numIndices, indices,
00245                            values, format, SUMINTO));
00246 }
00247 
00248 //----------------------------------------------------------------------------
00249 int Epetra_FECrsMatrix::SumIntoGlobalValues(int numRows, const int* rows,
00250               int numCols, const int* cols,
00251               const double* const* values,
00252               int format)
00253 {
00254   return(InputGlobalValues(numRows, rows,
00255                            numCols, cols,
00256                            values, format, SUMINTO));
00257 }
00258 
00259 //----------------------------------------------------------------------------
00260 int Epetra_FECrsMatrix::SumIntoGlobalValues(int numIndices, const int* indices,
00261               const double* values,
00262               int format)
00263 {
00264   return(InputGlobalValues(numIndices, indices,
00265                            numIndices, indices,
00266                            values, format, SUMINTO));
00267 }
00268 
00269 //----------------------------------------------------------------------------
00270 int Epetra_FECrsMatrix::SumIntoGlobalValues(int numRows, const int* rows,
00271               int numCols, const int* cols,
00272               const double* values,
00273               int format)
00274 {
00275   return(InputGlobalValues(numRows, rows,
00276                            numCols, cols,
00277                            values, format, SUMINTO));
00278 }
00279 
00280 //----------------------------------------------------------------------------
00281 int Epetra_FECrsMatrix::SumIntoGlobalValues(const Epetra_IntSerialDenseVector& indices,
00282               const Epetra_SerialDenseMatrix& values,
00283               int format)
00284 {
00285   if (indices.Length() != values.M() || indices.Length() != values.N()) {
00286     return(-1);
00287   }
00288 
00289   return( SumIntoGlobalValues(indices.Length(), indices.Values(),
00290             values.A(), format) );
00291 }
00292 
00293 //----------------------------------------------------------------------------
00294 int Epetra_FECrsMatrix::InsertGlobalValues(const Epetra_IntSerialDenseVector& indices,
00295               const Epetra_SerialDenseMatrix& values,
00296               int format)
00297 {
00298   if (indices.Length() != values.M() || indices.Length() != values.N()) {
00299     return(-1);
00300   }
00301 
00302   return( InsertGlobalValues(indices.Length(), indices.Values(),
00303             values.A(), format) );
00304 }
00305 
00306 //----------------------------------------------------------------------------
00307 int Epetra_FECrsMatrix::ReplaceGlobalValues(const Epetra_IntSerialDenseVector& indices,
00308               const Epetra_SerialDenseMatrix& values,
00309               int format)
00310 {
00311   if (indices.Length() != values.M() || indices.Length() != values.N()) {
00312     return(-1);
00313   }
00314 
00315   return( ReplaceGlobalValues(indices.Length(), indices.Values(),
00316             values.A(), format) );
00317 }
00318 
00319 //----------------------------------------------------------------------------
00320 int Epetra_FECrsMatrix::SumIntoGlobalValues(const Epetra_IntSerialDenseVector& rows,
00321               const Epetra_IntSerialDenseVector& cols,
00322               const Epetra_SerialDenseMatrix& values,
00323               int format)
00324 {
00325   if (rows.Length() != values.M() || cols.Length() != values.N()) {
00326     return(-1);
00327   }
00328 
00329   return( SumIntoGlobalValues(rows.Length(), rows.Values(),
00330             cols.Length(), cols.Values(),
00331             values.A(), format) );
00332 }
00333 
00334 //----------------------------------------------------------------------------
00335 int Epetra_FECrsMatrix::InsertGlobalValues(const Epetra_IntSerialDenseVector& rows,
00336              const Epetra_IntSerialDenseVector& cols,
00337              const Epetra_SerialDenseMatrix& values,
00338              int format)
00339 {
00340   if (rows.Length() != values.M() || cols.Length() != values.N()) {
00341     return(-1);
00342   }
00343 
00344   return( InsertGlobalValues(rows.Length(), rows.Values(),
00345            cols.Length(), cols.Values(),
00346             values.A(), format) );
00347 }
00348 
00349 //----------------------------------------------------------------------------
00350 int Epetra_FECrsMatrix::ReplaceGlobalValues(const Epetra_IntSerialDenseVector& rows,
00351               const Epetra_IntSerialDenseVector& cols,
00352               const Epetra_SerialDenseMatrix& values,
00353               int format)
00354 {
00355   if (rows.Length() != values.M() || cols.Length() != values.N()) {
00356     return(-1);
00357   }
00358 
00359   return( ReplaceGlobalValues(rows.Length(), rows.Values(),
00360             cols.Length(), cols.Values(),
00361             values.A(), format) );
00362 }
00363 
00364 //----------------------------------------------------------------------------
00365 int Epetra_FECrsMatrix::InsertGlobalValues(int numIndices, const int* indices,
00366             const double* const* values,
00367             int format)
00368 {
00369   return(InputGlobalValues(numIndices, indices,
00370                            numIndices, indices,
00371                            values, format, INSERT));
00372 }
00373 
00374 //----------------------------------------------------------------------------
00375 int Epetra_FECrsMatrix::InsertGlobalValues(int numRows, const int* rows,
00376             int numCols, const int* cols,
00377             const double* const* values,
00378             int format)
00379 {
00380   return(InputGlobalValues(numRows, rows,
00381                            numCols, cols,
00382                            values, format, INSERT));
00383 }
00384 
00385 //----------------------------------------------------------------------------
00386 int Epetra_FECrsMatrix::InsertGlobalValues(int numIndices, const int* indices,
00387               const double* values,
00388               int format)
00389 {
00390   return(InputGlobalValues(numIndices, indices,
00391                            numIndices, indices,
00392                            values, format, INSERT));
00393 }
00394 
00395 //----------------------------------------------------------------------------
00396 int Epetra_FECrsMatrix::InsertGlobalValues(int numRows, const int* rows,
00397               int numCols, const int* cols,
00398               const double* values,
00399               int format)
00400 {
00401   return(InputGlobalValues(numRows, rows,
00402                            numCols, cols,
00403                            values, format, INSERT));
00404 }
00405 
00406 //----------------------------------------------------------------------------
00407 int Epetra_FECrsMatrix::SumIntoGlobalValues(int GlobalRow, int NumEntries,
00408                                             const double* values, const int* Indices)
00409 {
00410   if (Map().MyGID(GlobalRow))
00411     return Epetra_CrsMatrix::SumIntoGlobalValues(GlobalRow, NumEntries,
00412             values, Indices);
00413   else if (useNonlocalMatrix_)
00414     return nonlocalMatrix_->SumIntoGlobalValues(GlobalRow,
00415            NumEntries, values, Indices);
00416   else
00417     return InputNonlocalGlobalValues(GlobalRow, NumEntries, Indices, values, SUMINTO);
00418 }
00419 
00420 //----------------------------------------------------------------------------
00421 int Epetra_FECrsMatrix::InsertGlobalValues(int GlobalRow, int NumEntries,
00422                                             const double* values, const int* Indices)
00423 {
00424   return(InputGlobalValues(1, &GlobalRow,
00425                            NumEntries, Indices, values,
00426                            Epetra_FECrsMatrix::ROW_MAJOR, INSERT));
00427 }
00428 
00429 //----------------------------------------------------------------------------
00430 int Epetra_FECrsMatrix::InsertGlobalValues(int GlobalRow, int NumEntries,
00431                                             double* values, int* Indices)
00432 {
00433   return(InputGlobalValues(1, &GlobalRow,
00434                            NumEntries, Indices, values,
00435                            Epetra_FECrsMatrix::ROW_MAJOR, INSERT));
00436 }
00437 
00438 //----------------------------------------------------------------------------
00439 int Epetra_FECrsMatrix::ReplaceGlobalValues(int GlobalRow, int NumEntries,
00440                                             const double* values, const int* Indices)
00441 {
00442   return(InputGlobalValues(1, &GlobalRow,
00443                            NumEntries, Indices, values,
00444                            Epetra_FECrsMatrix::ROW_MAJOR, REPLACE));
00445 }
00446 
00447 //----------------------------------------------------------------------------
00448 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numIndices, const int* indices,
00449               const double* const* values,
00450               int format)
00451 {
00452   return(InputGlobalValues(numIndices, indices,
00453                            numIndices, indices,
00454                            values, format, REPLACE));
00455 }
00456 
00457 //----------------------------------------------------------------------------
00458 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numRows, const int* rows,
00459               int numCols, const int* cols,
00460               const double* const* values,
00461               int format)
00462 {
00463   return(InputGlobalValues(numRows, rows,
00464                            numCols, cols,
00465                            values, format, REPLACE));
00466 }
00467 
00468 //----------------------------------------------------------------------------
00469 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numIndices, const int* indices,
00470               const double* values,
00471               int format)
00472 {
00473   return(InputGlobalValues(numIndices, indices,
00474                            numIndices, indices,
00475                            values, format, REPLACE));
00476 }
00477 
00478 //----------------------------------------------------------------------------
00479 int Epetra_FECrsMatrix::ReplaceGlobalValues(int numRows, const int* rows,
00480               int numCols, const int* cols,
00481               const double* values,
00482               int format)
00483 {
00484   return(InputGlobalValues(numRows, rows,
00485                            numCols, cols,
00486                            values, format, REPLACE));
00487 }
00488 
00489 //----------------------------------------------------------------------------
00490 int Epetra_FECrsMatrix::GlobalAssemble(bool callFillComplete, Epetra_CombineMode combineMode)
00491 {
00492   return( GlobalAssemble(DomainMap(), RangeMap(), callFillComplete, combineMode ));
00493 }
00494 
00495 //----------------------------------------------------------------------------
00496 int Epetra_FECrsMatrix::GlobalAssemble(const Epetra_Map& domain_map,
00497                                        const Epetra_Map& range_map,
00498                                        bool callFillComplete,
00499                                        Epetra_CombineMode combineMode)
00500 {
00501   if (Map().Comm().NumProc() < 2 || ignoreNonLocalEntries_) {
00502     if (callFillComplete) {
00503       EPETRA_CHK_ERR( FillComplete(domain_map, range_map) );
00504     }
00505     return(0);
00506   }
00507 
00508   Epetra_CrsMatrix* tempMat;
00509   if (useNonlocalMatrix_) {
00510     tempMat = nonlocalMatrix_;
00511   }
00512   else {
00513     //In this method we need to gather all the non-local (overlapping) data
00514     //that's been input on each processor, into the
00515     //non-overlapping distribution defined by the map that 'this' matrix was
00516     //constructed with.
00517 
00518     //First build a map that describes our nonlocal data.
00519     //We'll use the arbitrary distribution constructor of Map.
00520 
00521     int* nlr_ptr = nonlocalRows_.size() > 0 ? &nonlocalRows_[0] : 0;
00522     Epetra_Map* sourceMap = new Epetra_Map(-1, nonlocalRows_.size(), nlr_ptr,
00523             Map().IndexBase(), Map().Comm());
00524 
00525     //If sourceMap has global size 0, then no nonlocal data exists and we can
00526     //skip most of this function.
00527     if (sourceMap->NumGlobalElements() < 1) {
00528       if (callFillComplete) {
00529         EPETRA_CHK_ERR( FillComplete(domain_map, range_map) );
00530       }
00531       delete sourceMap;
00532       return(0);
00533     }
00534 
00535     //We also need to build a column-map, containing the columns in our
00536     //nonlocal data. To do that, create a list of all column-indices that
00537     //occur in our nonlocal rows.
00538 
00539     std::vector<int> cols;
00540 
00541     for(size_t i=0; i<nonlocalRows_.size(); ++i) {
00542       for(size_t j=0; j<nonlocalCols_[i].size(); ++j) {
00543         int col = nonlocalCols_[i][j];
00544         std::vector<int>::iterator it =
00545            std::lower_bound(cols.begin(), cols.end(), col);
00546         if (it == cols.end() || *it != col) {
00547           cols.insert(it, col);
00548         }
00549       }
00550     }
00551 
00552     int* cols_ptr = cols.size() > 0 ? &cols[0] : 0;
00553     Epetra_Map* colMap = new Epetra_Map(-1, cols.size(), cols_ptr,
00554          Map().IndexBase(), Map().Comm());
00555 
00556     //now we need to create a matrix with sourceMap and colMap, and fill it with
00557     //our nonlocal data so we can then export it to the correct owning processors.
00558 
00559     std::vector<int> nonlocalRowLengths(nonlocalRows_.size());
00560     for(size_t i=0; i<nonlocalRows_.size(); ++i) {
00561       nonlocalRowLengths[i] = nonlocalCols_[i].size();
00562     }
00563 
00564     int* nlRLptr = nonlocalRowLengths.size()>0 ? &nonlocalRowLengths[0] : NULL;
00565     tempMat = new Epetra_CrsMatrix(Copy, *sourceMap, *colMap, nlRLptr);
00566 
00567 
00568     for(size_t i=0; i<nonlocalRows_.size(); ++i) {
00569       EPETRA_CHK_ERR( tempMat->InsertGlobalValues(nonlocalRows_[i],
00570              nonlocalCols_[i].size(),
00571              &nonlocalCoefs_[i][0],
00572              &nonlocalCols_[i][0]) );
00573     }
00574 
00575     delete sourceMap;
00576     delete colMap;
00577   }
00578 
00579   //Next we need to make sure the 'indices-are-global' attribute of tempMat's
00580   //graph is set to true, in case this processor doesn't end up calling the
00581   //InsertGlobalValues method...
00582 
00583   const Epetra_CrsGraph& graph = tempMat->Graph();
00584   Epetra_CrsGraph& nonconst_graph = const_cast<Epetra_CrsGraph&>(graph);
00585   nonconst_graph.SetIndicesAreGlobal(true);
00586 
00587   //Now we need to call FillComplete on our temp matrix. We need to
00588   //pass a DomainMap and RangeMap, which are not the same as the RowMap
00589   //and ColMap that we constructed the matrix with.
00590   EPETRA_CHK_ERR(tempMat->FillComplete(domain_map, range_map));
00591 
00592   Epetra_Export* exporter = new Epetra_Export(tempMat->RowMap(), RowMap());
00593 
00594   EPETRA_CHK_ERR(Export(*tempMat, *exporter, combineMode));
00595 
00596   if(callFillComplete) {
00597     EPETRA_CHK_ERR(FillComplete(domain_map, range_map));
00598   }
00599 
00600   //now reset the values in our nonlocal data
00601   if (!useNonlocalMatrix_) {
00602     for(size_t i=0; i<nonlocalRows_.size(); ++i) {
00603       nonlocalCols_[i].resize(0);
00604       nonlocalCoefs_[i].resize(0);
00605     }
00606 
00607     delete tempMat;
00608   }
00609 
00610   delete exporter;
00611 
00612   return(0);
00613 }
00614 
00615 //----------------------------------------------------------------------------
00616 int Epetra_FECrsMatrix::InputGlobalValues_RowMajor(
00617             int numRows, const int* rows,
00618             int numCols, const int* cols,
00619             const double* values,
00620             int mode)
00621 {
00622   int returncode = 0;
00623   int err = 0;
00624 
00625   for(int i=0; i<numRows; ++i) {
00626     double* valuesptr = (double*)values + i*numCols;
00627 
00628     int local_row_id = Map().LID(rows[i]);
00629     if (local_row_id >= 0) {
00630       switch(mode) {
00631         case Epetra_FECrsMatrix::SUMINTO:
00632           err = this->Epetra_CrsMatrix::SumIntoGlobalValues(rows[i], numCols,
00633               valuesptr, (int*)cols);
00634           if (err<0) return(err);
00635           if (err>0) returncode = err;
00636           break;
00637         case Epetra_FECrsMatrix::REPLACE:
00638           err = this->Epetra_CrsMatrix::ReplaceGlobalValues(rows[i], numCols,
00639               valuesptr, (int*)cols);
00640           if (err<0) return(err);
00641           if (err>0) returncode = err;
00642           break;
00643         case Epetra_FECrsMatrix::INSERT:
00644           err = this->Epetra_CrsMatrix::InsertGlobalValues(rows[i], numCols,
00645               valuesptr, (int*)cols);
00646           if (err<0) return(err);
00647           if (err>0) returncode = err;
00648           break;
00649         default:
00650           cerr << "Epetra_FECrsMatrix: internal error, bad input mode."<<endl;
00651           return(-1);
00652       }
00653     }
00654     else {
00655       err = InputNonlocalGlobalValues(rows[i], numCols, cols,
00656           valuesptr, mode);
00657       if (err<0) return(err);
00658       if (err>0) returncode = err;
00659     }
00660   }
00661 
00662   return(returncode);
00663 }
00664 
00665 //----------------------------------------------------------------------------
00666 int Epetra_FECrsMatrix::InputGlobalValues(int numRows, const int* rows,
00667             int numCols, const int* cols,
00668             const double*const* values,
00669             int format, int mode)
00670 {
00671   if (format != Epetra_FECrsMatrix::ROW_MAJOR &&
00672       format != Epetra_FECrsMatrix::COLUMN_MAJOR) {
00673     cerr << "Epetra_FECrsMatrix: unrecognized format specifier."<<endl;
00674     return(-1);
00675   }
00676 
00677   if (format == Epetra_FECrsMatrix::COLUMN_MAJOR) {
00678     workData_.resize(numCols);
00679   }
00680 
00681   int returncode = 0;
00682 
00683   for(int i=0; i<numRows; ++i) {
00684     if (format == Epetra_FECrsMatrix::ROW_MAJOR) {
00685       returncode += InputGlobalValues_RowMajor(1, &rows[i], numCols, cols,
00686                                           values[i], mode);
00687       if (returncode < 0) return returncode;
00688       continue;
00689     }
00690 
00691     //If we get to here, the data is in column-major order.
00692 
00693     double* valuesptr = &workData_[0];
00694 
00695     //Since the data is in column-major order, then we copy the i-th row
00696     //of the values table into workData_, in order to have the row in
00697     //contiguous memory.
00698     //This is slow and not thread-safe.
00699 
00700     for(int j=0; j<numCols; ++j) {
00701       valuesptr[j] = values[j][i];
00702     }
00703 
00704     returncode += InputGlobalValues_RowMajor(1, &rows[i], numCols, cols, valuesptr, mode);
00705     if (returncode < 0) return returncode;
00706   }
00707 
00708   return(returncode);
00709 }
00710 
00711 //----------------------------------------------------------------------------
00712 int Epetra_FECrsMatrix::InputGlobalValues(int numRows, const int* rows,
00713             int numCols, const int* cols,
00714             const double* values,
00715             int format, int mode)
00716 {
00717   if (format == Epetra_FECrsMatrix::ROW_MAJOR) {
00718     return InputGlobalValues_RowMajor(numRows, rows, numCols, cols, values, mode);
00719   }
00720 
00721   workData_.resize(numCols);
00722 
00723   int returncode = 0;
00724   for(int i=0; i<numRows; ++i) {
00725     //copy each row out of the column-major values array, so we can pass it
00726     //to a row-major input function.
00727     for(int j=0; j<numCols; ++j) {
00728       workData_[j] = values[i+j*numRows];
00729     }
00730     int err = InputGlobalValues_RowMajor(1, &rows[i], numCols, cols, &workData_[0], mode);
00731     if (err < 0) return err;
00732     returncode += err;
00733   }
00734 
00735   return(returncode);
00736 }
00737 
00738 //----------------------------------------------------------------------------
00739 int Epetra_FECrsMatrix::InputNonlocalGlobalValues(int row,
00740               int numCols, const int* cols,
00741               const double* values,
00742               int mode)
00743 {
00744   // if we already have a nonlocal matrix object, this is easier...
00745   if (useNonlocalMatrix_) {
00746     int err, returncode = 0;
00747     double* valuesptr = (double*)values;
00748     switch(mode) {
00749     case Epetra_FECrsMatrix::SUMINTO:
00750       err = nonlocalMatrix_->SumIntoGlobalValues(row, numCols,
00751             valuesptr, (int*)cols);
00752       if (err<0) return(err);
00753       if (err>0) returncode = err;
00754       break;
00755     case Epetra_FECrsMatrix::REPLACE:
00756       err = nonlocalMatrix_->ReplaceGlobalValues(row, numCols,
00757             valuesptr, (int*)cols);
00758       if (err<0) return(err);
00759       if (err>0) returncode = err;
00760       break;
00761     case Epetra_FECrsMatrix::INSERT:
00762       err = nonlocalMatrix_->InsertGlobalValues(row, numCols,
00763            valuesptr, (int*)cols);
00764       if (err<0) return(err);
00765       if (err>0) returncode = err;
00766       break;
00767     default:
00768       cerr << "Epetra_FECrsMatrix: internal error, bad input mode."<<endl;
00769       return(-1);
00770     }
00771     return (returncode);
00772   }
00773 
00774   //find offset of this row in our list of nonlocal rows.
00775   std::vector<int>::iterator it =
00776       std::lower_bound(nonlocalRows_.begin(), nonlocalRows_.end(), row);
00777 
00778   int rowoffset = it - nonlocalRows_.begin();
00779   if (it == nonlocalRows_.end() || *it != row) {
00780     EPETRA_CHK_ERR( InsertNonlocalRow(row, it) );
00781   }
00782 
00783   for(int i=0; i<numCols; ++i) {
00784     EPETRA_CHK_ERR( InputNonlocalValue(rowoffset, cols[i], values[i], mode) );
00785   }
00786 
00787   return(0);
00788 }
00789 
00790 //----------------------------------------------------------------------------
00791 int Epetra_FECrsMatrix::InsertNonlocalRow(int row, std::vector<int>::iterator iter)
00792 {
00793   int offset = iter - nonlocalRows_.begin();
00794   nonlocalRows_.insert(iter, row);
00795   std::vector<std::vector<int> >::iterator cols_iter = nonlocalCols_.begin() + offset;
00796   nonlocalCols_.insert(cols_iter, std::vector<int>());
00797   std::vector<std::vector<double> >::iterator coefs_iter = nonlocalCoefs_.begin() + offset;
00798   nonlocalCoefs_.insert(coefs_iter, std::vector<double>());
00799 
00800   return(0);
00801 }
00802 
00803 //----------------------------------------------------------------------------
00804 int Epetra_FECrsMatrix::InputNonlocalValue(int rowoffset,
00805              int col, double value,
00806              int mode)
00807 {
00808   std::vector<int>& colIndices = nonlocalCols_[rowoffset];
00809   std::vector<double>& coefs = nonlocalCoefs_[rowoffset];
00810 
00811   std::vector<int>::iterator it =
00812      std::lower_bound(colIndices.begin(), colIndices.end(), col);
00813 
00814   if (it == colIndices.end() || *it != col) {
00815     int offset = it - colIndices.begin();
00816     colIndices.insert(it, col);
00817     std::vector<double>::iterator dit = coefs.begin()+offset;
00818     coefs.insert(dit, value);
00819     return 0;
00820   }
00821 
00822   int coloffset = it - colIndices.begin();
00823   if (mode == SUMINTO || mode == INSERT) {
00824     coefs[coloffset] += value;
00825   }
00826   else {
00827     coefs[coloffset] = value;
00828   }
00829 
00830   return(0);
00831 }
00832 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines