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