Epetra Package Browser (Single Doxygen Collection) Development
Epetra_VbrMatrix.h
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 
00044 #ifndef EPETRA_VBRMATRIX_H
00045 #define EPETRA_VBRMATRIX_H
00046 
00047 #include <Epetra_ConfigDefs.h>
00048 #include <Epetra_DistObject.h> 
00049 #include <Epetra_CompObject.h> 
00050 #include <Epetra_BLAS.h>
00051 #include <Epetra_RowMatrix.h>
00052 #include <Epetra_Operator.h>
00053 #include <Epetra_CrsGraph.h>
00054 #include <Epetra_SerialDenseMatrix.h>
00055 class Epetra_BlockMap;
00056 class Epetra_Map;
00057 class Epetra_Import;
00058 class Epetra_Export;
00059 class Epetra_Vector;
00060 class Epetra_MultiVector;
00061 
00063 
00170 class EPETRA_LIB_DLL_EXPORT Epetra_VbrMatrix : public Epetra_DistObject,
00171           public Epetra_CompObject,
00172           public Epetra_BLAS,
00173           public virtual Epetra_RowMatrix {
00174  public:
00175 
00177 
00178 
00179 
00189   Epetra_VbrMatrix(Epetra_DataAccess CV, const Epetra_BlockMap& RowMap, int *NumBlockEntriesPerRow);
00190   
00192 
00203   Epetra_VbrMatrix(Epetra_DataAccess CV, const Epetra_BlockMap& RowMap, int NumBlockEntriesPerRow);
00204 
00206 
00218   Epetra_VbrMatrix(Epetra_DataAccess CV, const Epetra_BlockMap& RowMap, const Epetra_BlockMap& ColMap, int *NumBlockEntriesPerRow);
00219   
00221 
00234   Epetra_VbrMatrix(Epetra_DataAccess CV, const Epetra_BlockMap& RowMap, const Epetra_BlockMap& ColMap, int NumBlockEntriesPerRow);
00235 
00237 
00246   Epetra_VbrMatrix(Epetra_DataAccess CV, const Epetra_CrsGraph & Graph);
00247 
00249   Epetra_VbrMatrix(const Epetra_VbrMatrix & Matrix);
00250 
00252   virtual ~Epetra_VbrMatrix();
00254   
00256 
00257 
00258   Epetra_VbrMatrix& operator=(const Epetra_VbrMatrix& src);
00259 
00261 
00267     int PutScalar(double ScalarConstant);
00268 
00270 
00276     int Scale(double ScalarConstant);
00277 
00278 
00280 
00290     int BeginInsertGlobalValues(int BlockRow,
00291         int NumBlockEntries,
00292         int * BlockIndices);
00293 
00295 
00305     int BeginInsertMyValues(int BlockRow, int NumBlockEntries, int * BlockIndices);
00306 
00308 
00318     int BeginReplaceGlobalValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00319 
00321 
00331     int BeginReplaceMyValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00332 
00334 
00344     int BeginSumIntoGlobalValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00345 
00347 
00357     int BeginSumIntoMyValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00358 
00360     /* Submit a block entry that will recorded in the block row that was initiated by one of the
00361        Begin routines listed above.  Once a one of the following routines: BeginInsertGlobalValues(),
00362        BeginInsertMyValues(), BeginReplaceGlobalValues(), BeginReplaceMyValues(), BeginSumIntoGlobalValues(),
00363        BeginSumIntoMyValues(), you \e must call SubmitBlockEntry() NumBlockEntries times to register the values 
00364        corresponding to the block indices passed in to the Begin routine.  If the Epetra_VbrMatrix constuctor
00365        was called in Copy mode, the values will be copied.  However, no copying will be done until the EndSubmitEntries()
00366        function is call to complete submission of the current block row.  If the constructor was called in View mode, all
00367        block entries passed via SubmitBlockEntry() will not be copied, but a pointer will be set to point to the argument Values
00368        that was passed in by the user.
00369 
00370        For performance reasons, SubmitBlockEntry() does minimal processing of data.  Any processing that can be
00371        delayed is performed in EndSubmitEntries().
00372 
00373     \param In
00374            Values - The starting address of the values.
00375     \param In
00376            LDA - The stride between successive columns of Values.
00377     \param In
00378            NumRows - The number of rows passed in.
00379     \param In
00380            NumCols - The number of columns passed in.
00381 
00382     \return Integer error code, set to 0 if successful.
00383     */
00384     int SubmitBlockEntry(double *Values, int LDA, int NumRows, int NumCols);
00385 
00387     /* Submit a block entry that will recorded in the block row that was initiated by one of the
00388        Begin routines listed above.  Once a one of the following routines: BeginInsertGlobalValues(),
00389        BeginInsertMyValues(), BeginReplaceGlobalValues(), BeginReplaceMyValues(), BeginSumIntoGlobalValues(),
00390        BeginSumIntoMyValues(), you \e must call SubmitBlockEntry() NumBlockEntries times to register the values 
00391        corresponding to the block indices passed in to the Begin routine.  If the Epetra_VbrMatrix constuctor
00392        was called in Copy mode, the values will be copied.  However, no copying will be done until the EndSubmitEntries()
00393        function is call to complete submission of the current block row.  If the constructor was called in View mode, all
00394        block entries passed via SubmitBlockEntry() will not be copied, but a pointer will be set to point to the argument Values
00395        that was passed in by the user.
00396 
00397        For performance reasons, SubmitBlockEntry() does minimal processing of data.  Any processing that can be
00398        delayed is performed in EndSubmitEntries().
00399 
00400     \param In
00401            Mat - Preformed dense matrix block.
00402 
00403     \return Integer error code, set to 0 if successful.
00404     */
00405     int SubmitBlockEntry( Epetra_SerialDenseMatrix &Mat );
00406 
00408 
00413     int EndSubmitEntries();
00414 
00416 
00427     int ReplaceDiagonalValues(const Epetra_Vector & Diagonal);
00428 
00430     /* This version of FillComplete assumes that the domain and range
00431        distributions are identical to the matrix row distributions.
00432        \return error code, 0 if successful. Returns a positive warning code of 3
00433         if the matrix is rectangular (meaning that the other overloading of
00434         FillComplete should have been called, with differen domain-map and
00435         range-map specified).
00436     */
00437     int FillComplete();
00438 
00440     /* This version of FillComplete requires the explicit specification of the domain
00441        and range distribution maps.  These maps are used for importing and exporting vector
00442        and multi-vector elements that are needed for distributed matrix computations.  For
00443        example, to compute y = Ax in parallel, we would specify the DomainMap as the distribution
00444        of the vector x and the RangeMap as the distribution of the vector y.
00445     \param In
00446            DomainMap - Map that describes the distribution of vector and multi-vectors in the
00447                  matrix domain.
00448     \param In
00449            RangeMap - Map that describes the distribution of vector and multi-vectors in the
00450                  matrix range.
00451 
00452     \return error code, 0 if successful. positive warning code of 2 if it is detected that the
00453     matrix-graph got out of sync since this matrix was constructed (for instance if
00454     graph.FillComplete() was called by another matrix that shares the graph)
00455     */
00456     int FillComplete(const Epetra_BlockMap& DomainMap, const Epetra_BlockMap& RangeMap);
00457 
00459     bool Filled() const {return(Graph_->Filled());};
00461 
00463 
00464 
00466 
00488     int ExtractGlobalBlockRowPointers(int BlockRow, int MaxNumBlockEntries, 
00489               int & RowDim,  int & NumBlockEntries, 
00490               int * BlockIndices,
00491               Epetra_SerialDenseMatrix ** & Values) const;
00492 
00494 
00516     int ExtractMyBlockRowPointers(int BlockRow, int MaxNumBlockEntries, 
00517           int & RowDim, int & NumBlockEntries, 
00518           int * BlockIndices,
00519           Epetra_SerialDenseMatrix** & Values) const;
00520 
00522 
00538     int BeginExtractGlobalBlockRowCopy(int BlockRow, int MaxNumBlockEntries, 
00539                int & RowDim,  int & NumBlockEntries, 
00540                int * BlockIndices, int * ColDims) const;
00541 
00543 
00559     int BeginExtractMyBlockRowCopy(int BlockRow, int MaxNumBlockEntries, 
00560            int & RowDim, int & NumBlockEntries, 
00561            int * BlockIndices, int * ColDims) const;
00562 
00564 
00579     int ExtractEntryCopy(int SizeOfValues, double * Values, int LDA, bool SumInto) const;
00580 
00582 
00594     int BeginExtractGlobalBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00595                int * & BlockIndices) const;
00596 
00598 
00610     int BeginExtractMyBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00611            int * & BlockIndices) const;
00612 
00613 
00615 
00622     int ExtractEntryView(Epetra_SerialDenseMatrix* & entry) const;
00623 
00625 
00639     int ExtractGlobalBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00640           int * & BlockIndices,
00641           Epetra_SerialDenseMatrix** & Values) const;
00642 
00644 
00658     int ExtractMyBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00659             int * & BlockIndices,
00660             Epetra_SerialDenseMatrix** & Values) const;
00661 
00662 
00664 
00670     int ExtractDiagonalCopy(Epetra_Vector & Diagonal) const;
00671 
00673 
00683     int BeginExtractBlockDiagonalCopy(int MaxNumBlockDiagonalEntries, 
00684               int & NumBlockDiagonalEntries, int * RowColDims ) const;
00686 
00701     int ExtractBlockDiagonalEntryCopy(int SizeOfValues, double * Values, int LDA, bool SumInto) const;
00702 
00704 
00712     int BeginExtractBlockDiagonalView(int & NumBlockDiagonalEntries, int * & RowColDims ) const;
00713 
00715 
00725     int ExtractBlockDiagonalEntryView(double * & Values, int & LDA) const;
00727 
00729 
00730 
00731 
00733 
00743     int Multiply1(bool TransA, const Epetra_Vector& x, Epetra_Vector& y) const;
00744 
00746 
00756     int Multiply(bool TransA, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00757 
00759 
00773     int Solve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_Vector& x, Epetra_Vector& y) const;
00774 
00776 
00790     int Solve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00791 
00792 
00794 
00804     int InvRowSums(Epetra_Vector& x) const;
00805 
00807 
00814     int LeftScale(const Epetra_Vector& x);
00815 
00817 
00827     int InvColSums(Epetra_Vector& x) const ;
00828 
00830 
00837     int RightScale(const Epetra_Vector& x);
00839 
00841 
00842 
00844     int OptimizeStorage();
00845 
00847     bool StorageOptimized() const {return(StorageOptimized_);};
00848 
00850     bool IndicesAreGlobal() const {return(Graph_->IndicesAreGlobal());};
00851 
00853     bool IndicesAreLocal() const {return(Graph_->IndicesAreLocal());};
00854 
00856     bool IndicesAreContiguous() const {return(Graph_->IndicesAreContiguous());};
00857 
00859     bool LowerTriangular() const {return(Graph_->LowerTriangular());};
00860 
00862     bool UpperTriangular() const {return(Graph_->UpperTriangular());};
00863 
00865     bool NoDiagonal() const {return(Graph_->NoDiagonal());};
00866 
00868   
00870 
00871 
00873     /* Returns the quantity \f$ \| A \|_\infty\f$ such that
00874        \f[\| A \|_\infty = \max_{1\lei\lem} \sum_{j=1}^n |a_{ij}| \f].
00875      \warning The NormInf() method will not properly calculate the infinity norm for a matrix that has entries that are
00876      replicated on multiple processors. */ 
00877     double NormInf() const;
00878 
00880     /* Returns the quantity \f$ \| A \|_1\f$ such that
00881        \f[\| A \|_1 = \max_{1\lej\len} \sum_{i=1}^m |a_{ij}| \f].
00882      \warning The NormOne() method will not properly calculate the one norm for a matrix that has entries that are
00883     */ 
00884     double NormOne() const;
00885 
00887     /* Returns the quantity \f[ \| A \|_{Frobenius} = \sqrt{\sum_{i=1}^m \sum_{j=1}^n\|a_{ij}\|^2}\f]
00888        \warning the NormFrobenius() method will not properly calculate the frobenius norm for a matrix that
00889        has entries which are replicated on multiple processors. In that case, the returned
00890        norm will be larger than the true norm.
00891     */
00892     double NormFrobenius() const;
00893 
00895     int MaxRowDim() const {return(Graph_->MaxRowDim());};
00896 
00898     int MaxColDim() const {return(Graph_->MaxColDim());};
00899 
00901     int GlobalMaxRowDim() const {return(Graph_->GlobalMaxRowDim());};
00902 
00904     int GlobalMaxColDim() const {return(Graph_->GlobalMaxColDim());};
00905 
00907     int NumMyRows() const {return(Graph_->NumMyRows());};
00909     int NumMyCols() const {return(Graph_->NumMyCols());};
00910 
00912     int NumMyNonzeros() const {return(Graph_->NumMyNonzeros());};
00913 
00915     int NumGlobalRows() const {return(Graph_->NumGlobalRows());};
00916 
00918     int NumGlobalCols() const {return(Graph_->NumGlobalCols());};
00919 
00921     /*
00922      Note that if maps are defined such that some nonzeros appear on
00923      multiple processors, then those nonzeros will be counted multiple
00924      times. If the user wishes to assemble a matrix from overlapping
00925      submatrices, they can use Epetra_FEVbrMatrix.
00926     */
00927     int NumGlobalNonzeros() const {return(Graph_->NumGlobalNonzeros());};
00928 
00930     int NumMyBlockRows() const {return(Graph_->NumMyBlockRows());};
00931 
00933     int NumMyBlockCols() const {return(Graph_->NumMyBlockCols());};
00934     
00936     int NumMyBlockEntries() const {return(Graph_->NumMyEntries());};
00937 
00939     int NumMyBlockDiagonals() const {return(Graph_->NumMyBlockDiagonals());};
00940     
00942     int NumMyDiagonals() const {return(Graph_->NumMyDiagonals());};
00943     
00945     int NumGlobalBlockRows() const {return(Graph_->NumGlobalBlockRows());};
00946     
00948     int NumGlobalBlockCols() const {return(Graph_->NumGlobalBlockCols());};
00949     
00951     int NumGlobalBlockEntries() const {return(Graph_->NumGlobalEntries());};
00952     
00954     int NumGlobalBlockDiagonals() const {return(Graph_->NumGlobalBlockDiagonals());};
00955     
00957     int NumGlobalDiagonals() const {return(Graph_->NumGlobalDiagonals());};
00958 
00960     int NumGlobalBlockEntries(int Row) const {return(Graph_->NumGlobalIndices(Row));};
00961 
00963     int NumAllocatedGlobalBlockEntries(int Row) const{return(Graph_->NumAllocatedGlobalIndices(Row));};
00964 
00966     int MaxNumBlockEntries() const {return(Graph_->MaxNumIndices());};
00967 
00969     int GlobalMaxNumBlockEntries() const {return(Graph_->GlobalMaxNumIndices());};
00970 
00972     int NumMyBlockEntries(int Row) const {return(Graph_->NumMyIndices(Row));};
00973 
00975     int NumAllocatedMyBlockEntries(int Row) const {return(Graph_->NumAllocatedMyIndices(Row));};
00976 
00978 
00982     int MaxNumNonzeros() const {return(Graph_->MaxNumNonzeros());};
00983 
00985 
00987     int GlobalMaxNumNonzeros() const {return(Graph_->GlobalMaxNumNonzeros());};
00988 
00990     int IndexBase() const {return(Graph_->IndexBase());};
00991 
00993     const Epetra_CrsGraph & Graph() const {return(*Graph_);};
00994 
00996     const Epetra_Import * Importer() const {return(Graph_->Importer());};
00997 
00999     const Epetra_Export * Exporter() const {return(Graph_->Exporter());};
01000 
01002     const Epetra_BlockMap & DomainMap() const {return(Graph_->DomainMap());};
01003 
01005     const Epetra_BlockMap & RangeMap() const  {return(Graph_->RangeMap());};
01006 
01008     const Epetra_BlockMap & RowMap() const {return(Graph_->RowMap());};
01009 
01011     const Epetra_BlockMap & ColMap() const {return(Graph_->ColMap());};
01012 
01014 
01016     const Epetra_Comm & Comm() const {return(Graph_->Comm());};
01017 
01019   
01021 
01022 
01023     int LRID( int GRID_in) const {return(Graph_->LRID(GRID_in));};
01024 
01026     int GRID( int LRID_in) const {return(Graph_->GRID(LRID_in));};
01027 
01029     int LCID( int GCID_in) const {return(Graph_->LCID(GCID_in));};
01030 
01032     int GCID( int LCID_in) const {return(Graph_->GCID(LCID_in));};
01033  
01035     bool  MyGRID(int GRID_in) const {return(Graph_->MyGRID(GRID_in));};
01036    
01038     bool  MyLRID(int LRID_in) const {return(Graph_->MyLRID(LRID_in));};
01039 
01041     bool  MyGCID(int GCID_in) const {return(Graph_->MyGCID(GCID_in));};
01042    
01044     bool  MyLCID(int LCID_in) const {return(Graph_->MyLCID(LCID_in));};
01045 
01047     bool MyGlobalBlockRow(int GID) const {return(Graph_->MyGlobalRow(GID));};
01049   
01051 
01052 
01054   virtual void Print(ostream & os) const;
01056 
01058 
01059 
01061     const char * Label() const {return(Epetra_Object::Label());};
01062     
01064 
01073   int SetUseTranspose(bool UseTranspose_in) {UseTranspose_ = UseTranspose_in; return(0);};
01074 
01076 
01084   int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
01085 
01087 
01100   int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
01101 
01103     bool HasNormInf() const {return(true);};
01104 
01106     bool UseTranspose() const {return(UseTranspose_);};
01107 
01109     const Epetra_Map & OperatorDomainMap() const 
01110     {
01111       if (!HavePointObjects_) GeneratePointObjects();
01112       if (UseTranspose()) return(*OperatorRangeMap_);
01113       else return(*OperatorDomainMap_);
01114     }
01115 
01117     const Epetra_Map & OperatorRangeMap() const 
01118     {
01119       if (!HavePointObjects_) GeneratePointObjects();
01120       if (UseTranspose()) return(*OperatorDomainMap_);
01121       else return(*OperatorRangeMap_);
01122     }
01123 
01125 
01126 
01127 
01129 
01143     int ExtractGlobalRowCopy(int GlobalRow, int Length, int & NumEntries, double *Values, int * Indices) const;
01144 
01146 
01160     int ExtractMyRowCopy(int MyRow, int Length, int & NumEntries, double *Values, int * Indices) const;
01161 
01163 
01171     int NumMyRowEntries(int MyRow, int & NumEntries) const;
01172 
01174     int MaxNumEntries() const;
01175 
01177     const Epetra_BlockMap& Map() const { return Epetra_DistObject::Map(); }
01178 
01180     const Epetra_Map & RowMatrixRowMap() const 
01181       { if (!HavePointObjects_) GeneratePointObjects(); return(*RowMatrixRowMap_); };
01182 
01184     const Epetra_Map & RowMatrixColMap() const 
01185       { if (!HavePointObjects_) GeneratePointObjects(); return(*RowMatrixColMap_); };
01186 
01188     const Epetra_Import * RowMatrixImporter() const 
01189       { if (!HavePointObjects_) GeneratePointObjects(); return(RowMatrixImporter_); };
01190 
01192 
01194 
01195 
01197     const Epetra_BlockMap & BlockImportMap() const {return(Graph_->ImportMap());};
01198 
01200     int TransformToLocal();
01201     
01203     int TransformToLocal(const Epetra_BlockMap* DomainMap, const Epetra_BlockMap* RangeMap);
01205 
01206  protected:
01207     void DeleteMemory();
01208     bool Allocated() const {return(Allocated_);};
01209     int SetAllocated(bool Flag) {Allocated_ = Flag; return(0);};
01210     Epetra_SerialDenseMatrix *** Values() const {return(Entries_);};
01211 
01212   // Internal utilities
01213 
01214   int DoMultiply(bool TransA, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
01215   int DoSolve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
01216   void InitializeDefaults();
01217   int Allocate();
01218   int BeginInsertValues(int BlockRow, int NumBlockEntries, 
01219       int * BlockIndices, bool IndicesAreLocal);
01220   int BeginReplaceValues(int BlockRow, int NumBlockEntries, 
01221        int *BlockIndices, bool IndicesAreLocal);
01222   int BeginSumIntoValues(int BlockRow, int NumBlockEntries, 
01223        int *BlockIndices, bool IndicesAreLocal);
01224   int SetupForSubmits(int BlockRow, int NumBlockEntries, int * BlockIndices, 
01225           bool IndicesAreLocal, Epetra_CombineMode SubmitMode);
01226   int EndReplaceSumIntoValues();
01227   int EndInsertValues();
01228 
01229   int CopyMat(double * A, int LDA, int NumRows, int NumCols, 
01230         double * B, int LDB, bool SumInto) const;
01231   int BeginExtractBlockRowCopy(int BlockRow, int MaxNumBlockEntries, 
01232              int & RowDim, int & NumBlockEntries, 
01233              int * BlockIndices, int * ColDims, 
01234              bool IndicesAreLocal) const;
01235   int SetupForExtracts(int BlockRow, int & RowDim, int NumBlockEntries,
01236            bool ExtractView, bool IndicesAreLocal) const;
01237   int ExtractBlockDimsCopy(int NumBlockEntries, int * ColDims) const;
01238   int ExtractBlockRowPointers(int BlockRow, int MaxNumBlockEntries, 
01239             int & RowDim, int & NumBlockEntries, 
01240             int * BlockIndices, 
01241             Epetra_SerialDenseMatrix ** & Values,
01242             bool IndicesAreLocal) const;
01243   int BeginExtractBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
01244              int * & BlockIndices, 
01245              bool IndicesAreLocal) const;
01246   int CopyMatDiag(double * A, int LDA, int NumRows, int NumCols, 
01247       double * Diagonal) const;
01248   int ReplaceMatDiag(double * A, int LDA, int NumRows, int NumCols, 
01249          double * Diagonal);
01250 
01251   //This BlockRowMultiply accepts Alpha and Beta arguments. It is called
01252   //from within the 'solve' methods.
01253   void BlockRowMultiply(bool TransA, int RowDim, int NumEntries, 
01254       int * BlockIndices, int RowOff,
01255       int * FirstPointInElementList, int * ElementSizeList,
01256       double Alpha, Epetra_SerialDenseMatrix** As,
01257       double ** X, double Beta, double ** Y, int NumVectors) const;
01258 
01259   //This BlockRowMultiply doesn't accept Alpha and Beta arguments, instead it
01260   //assumes that they are both 1.0. It is called from within the 'Multiply'
01261   //methods.
01262   void BlockRowMultiply(bool TransA, int RowDim, int NumEntries, 
01263       int * BlockIndices, int RowOff,
01264       int * FirstPointInElementList,
01265       int * ElementSizeList,
01266       Epetra_SerialDenseMatrix** As,
01267       double ** X, double ** Y, int NumVectors) const;
01268   //
01269   // Assumes Alpha=Beta=1 and works only on storage optimized matrices
01270   //
01271   void FastBlockRowMultiply(bool TransA, int RowDim, int NumEntries, 
01272           int * BlockIndices, int RowOff,
01273           int * FirstPointInElementList,
01274           int * ElementSizeList,
01275           Epetra_SerialDenseMatrix** As,
01276           double ** X, double ** Y, int NumVectors) const;
01277 
01278   int InverseSums(bool DoRows, Epetra_Vector& x) const;
01279   int Scale(bool DoRows, const Epetra_Vector& x);
01280   void BlockRowNormInf(int RowDim, int NumEntries, 
01281            Epetra_SerialDenseMatrix** As, 
01282            double * Y) const;
01283   void BlockRowNormOne(int RowDim, int NumEntries, int * BlockRowIndices,
01284            Epetra_SerialDenseMatrix** As, 
01285            int * ColFirstPointInElementList, double * x) const;
01286   void SetStaticGraph(bool Flag) {StaticGraph_ = Flag;};
01287 
01288   int CheckSizes(const Epetra_SrcDistObject& A);
01289 
01290   int CopyAndPermute(const Epetra_SrcDistObject & Source,
01291                      int NumSameIDs, 
01292                      int NumPermuteIDs,
01293                      int * PermuteToLIDs,
01294                      int *PermuteFromLIDs,
01295                      const Epetra_OffsetIndex * Indexor);
01296   
01297   int PackAndPrepare(const Epetra_SrcDistObject & Source,
01298                      int NumExportIDs,
01299                      int * ExportLIDs,
01300                      int & LenExports,
01301                      char * & Exports,
01302                      int & SizeOfPacket,
01303                      int * Sizes,
01304                      bool & VarSizes,
01305                      Epetra_Distributor & Distor);
01306   
01307   int UnpackAndCombine(const Epetra_SrcDistObject & Source, 
01308                        int NumImportIDs,
01309                        int * ImportLIDs, 
01310                        int LenImports,
01311                        char * Imports,
01312                        int & SizeOfPacket, 
01313                        Epetra_Distributor & Distor,
01314                        Epetra_CombineMode CombineMode,
01315                        const Epetra_OffsetIndex * Indexor);
01316 
01318   int SortEntries();
01319 
01321   bool Sorted() const {return(Graph_->Sorted());};
01322 
01324   int MergeRedundantEntries();
01325 
01327   bool NoRedundancies() const {return(Graph_->NoRedundancies());};
01328 
01329   bool StaticGraph() const {return(StaticGraph_);};
01330 
01331   int GeneratePointObjects() const;
01332   int BlockMap2PointMap(const Epetra_BlockMap & BlockMap, Epetra_Map * & PointMap) const;
01333   int UpdateOperatorXY(const Epetra_MultiVector& X, const Epetra_MultiVector& Y) const;
01334 
01335   Epetra_CrsGraph * Graph_;
01336   bool Allocated_;
01337   bool StaticGraph_;
01338   bool UseTranspose_;
01339   bool constructedWithFilledGraph_;
01340   bool matrixFillCompleteCalled_;
01341   bool StorageOptimized_;
01342 
01343   int NumMyBlockRows_;
01344 
01345   Epetra_DataAccess CV_;
01346 
01347 
01348   int * NumBlockEntriesPerRow_;
01349   int * NumAllocatedBlockEntriesPerRow_;
01350   int ** Indices_;
01351   int * ElementSizeList_;
01352   int * FirstPointInElementList_;
01353 
01354   Epetra_SerialDenseMatrix ***Entries_;
01355 
01356   double *All_Values_Orig_;
01357   double *All_Values_;
01358 
01359   mutable double NormInf_;
01360   mutable double NormOne_;
01361   mutable double NormFrob_;
01362 
01363   mutable Epetra_MultiVector * ImportVector_;
01364   mutable Epetra_MultiVector * ExportVector_;
01365 
01366   // State variables needed for constructing matrix entry-by-entry
01367   mutable int *TempRowDims_;
01368   mutable Epetra_SerialDenseMatrix **TempEntries_;
01369   mutable int LenTemps_;
01370   mutable int CurBlockRow_;
01371   mutable int CurNumBlockEntries_;
01372   mutable int * CurBlockIndices_;
01373   mutable int CurEntry_;
01374   mutable bool CurIndicesAreLocal_;
01375   mutable Epetra_CombineMode CurSubmitMode_;
01376   
01377   // State variables needed for extracting entries
01378   mutable int CurExtractBlockRow_;
01379   mutable int CurExtractEntry_; 
01380   mutable int CurExtractNumBlockEntries_;
01381   mutable bool CurExtractIndicesAreLocal_;
01382   mutable bool CurExtractView_;
01383   mutable int CurRowDim_;
01384 
01385   // State variable for extracting block diagonal entries
01386   mutable int CurBlockDiag_;
01387 
01388   // Maps and importer that support the Epetra_RowMatrix interface
01389   mutable Epetra_Map * RowMatrixRowMap_;
01390   mutable Epetra_Map * RowMatrixColMap_;
01391   mutable Epetra_Import * RowMatrixImporter_;
01392 
01393   // Maps that support the Epetra_Operator interface
01394   mutable Epetra_Map * OperatorDomainMap_;
01395   mutable Epetra_Map * OperatorRangeMap_;
01396   mutable Epetra_MultiVector * OperatorX_;
01397   mutable Epetra_MultiVector * OperatorY_;
01398 
01399   // bool to indicate if above four point maps and importer have already been created
01400   mutable bool HavePointObjects_;
01401 
01402   bool squareFillCompleteCalled_;
01403 };
01404 
01405 #endif /* EPETRA_VBRMATRIX_H */
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines