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 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 
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 
00279 
00284     int DirectSubmitBlockEntry(int GlobalBlockRow, int GlobalBlockCol,
00285                                const double *values, int LDA,
00286                                int NumRows, int NumCols, bool sum_into);
00287 
00289 
00299     int BeginInsertGlobalValues(int BlockRow,
00300         int NumBlockEntries,
00301         int * BlockIndices);
00302 
00304 
00314     int BeginInsertMyValues(int BlockRow, int NumBlockEntries, int * BlockIndices);
00315 
00317 
00327     int BeginReplaceGlobalValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00328 
00330 
00340     int BeginReplaceMyValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00341 
00343 
00353     int BeginSumIntoGlobalValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00354 
00356 
00366     int BeginSumIntoMyValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00367 
00369     /* Submit a block entry that will recorded in the block row that was initiated by one of the
00370        Begin routines listed above.  Once a one of the following routines: BeginInsertGlobalValues(),
00371        BeginInsertMyValues(), BeginReplaceGlobalValues(), BeginReplaceMyValues(), BeginSumIntoGlobalValues(),
00372        BeginSumIntoMyValues(), you \e must call SubmitBlockEntry() NumBlockEntries times to register the values 
00373        corresponding to the block indices passed in to the Begin routine.  If the Epetra_VbrMatrix constuctor
00374        was called in Copy mode, the values will be copied.  However, no copying will be done until the EndSubmitEntries()
00375        function is call to complete submission of the current block row.  If the constructor was called in View mode, all
00376        block entries passed via SubmitBlockEntry() will not be copied, but a pointer will be set to point to the argument Values
00377        that was passed in by the user.
00378 
00379        For performance reasons, SubmitBlockEntry() does minimal processing of data.  Any processing that can be
00380        delayed is performed in EndSubmitEntries().
00381 
00382     \param In
00383            Values - The starting address of the values.
00384     \param In
00385            LDA - The stride between successive columns of Values.
00386     \param In
00387            NumRows - The number of rows passed in.
00388     \param In
00389            NumCols - The number of columns passed in.
00390 
00391     \return Integer error code, set to 0 if successful.
00392     */
00393     int SubmitBlockEntry(double *Values, int LDA, int NumRows, int NumCols);
00394 
00396     /* Submit a block entry that will recorded in the block row that was initiated by one of the
00397        Begin routines listed above.  Once a one of the following routines: BeginInsertGlobalValues(),
00398        BeginInsertMyValues(), BeginReplaceGlobalValues(), BeginReplaceMyValues(), BeginSumIntoGlobalValues(),
00399        BeginSumIntoMyValues(), you \e must call SubmitBlockEntry() NumBlockEntries times to register the values 
00400        corresponding to the block indices passed in to the Begin routine.  If the Epetra_VbrMatrix constuctor
00401        was called in Copy mode, the values will be copied.  However, no copying will be done until the EndSubmitEntries()
00402        function is call to complete submission of the current block row.  If the constructor was called in View mode, all
00403        block entries passed via SubmitBlockEntry() will not be copied, but a pointer will be set to point to the argument Values
00404        that was passed in by the user.
00405 
00406        For performance reasons, SubmitBlockEntry() does minimal processing of data.  Any processing that can be
00407        delayed is performed in EndSubmitEntries().
00408 
00409     \param In
00410            Mat - Preformed dense matrix block.
00411 
00412     \return Integer error code, set to 0 if successful.
00413     */
00414     int SubmitBlockEntry( Epetra_SerialDenseMatrix &Mat );
00415 
00417 
00422     int EndSubmitEntries();
00423 
00425 
00436     int ReplaceDiagonalValues(const Epetra_Vector & Diagonal);
00437 
00439     /* This version of FillComplete assumes that the domain and range
00440        distributions are identical to the matrix row distributions.
00441        \return error code, 0 if successful. Returns a positive warning code of 3
00442         if the matrix is rectangular (meaning that the other overloading of
00443         FillComplete should have been called, with differen domain-map and
00444         range-map specified).
00445     */
00446     int FillComplete();
00447 
00449     /* This version of FillComplete requires the explicit specification of the domain
00450        and range distribution maps.  These maps are used for importing and exporting vector
00451        and multi-vector elements that are needed for distributed matrix computations.  For
00452        example, to compute y = Ax in parallel, we would specify the DomainMap as the distribution
00453        of the vector x and the RangeMap as the distribution of the vector y.
00454     \param In
00455            DomainMap - Map that describes the distribution of vector and multi-vectors in the
00456                  matrix domain.
00457     \param In
00458            RangeMap - Map that describes the distribution of vector and multi-vectors in the
00459                  matrix range.
00460 
00461     \return error code, 0 if successful. positive warning code of 2 if it is detected that the
00462     matrix-graph got out of sync since this matrix was constructed (for instance if
00463     graph.FillComplete() was called by another matrix that shares the graph)
00464     */
00465     int FillComplete(const Epetra_BlockMap& DomainMap, const Epetra_BlockMap& RangeMap);
00466 
00468     bool Filled() const {return(Graph_->Filled());};
00470 
00472 
00473 
00475 
00497     int ExtractGlobalBlockRowPointers(int BlockRow, int MaxNumBlockEntries, 
00498               int & RowDim,  int & NumBlockEntries, 
00499               int * BlockIndices,
00500               Epetra_SerialDenseMatrix ** & Values) const;
00501 
00503 
00525     int ExtractMyBlockRowPointers(int BlockRow, int MaxNumBlockEntries, 
00526           int & RowDim, int & NumBlockEntries, 
00527           int * BlockIndices,
00528           Epetra_SerialDenseMatrix** & Values) const;
00529 
00531 
00547     int BeginExtractGlobalBlockRowCopy(int BlockRow, int MaxNumBlockEntries, 
00548                int & RowDim,  int & NumBlockEntries, 
00549                int * BlockIndices, int * ColDims) const;
00550 
00552 
00568     int BeginExtractMyBlockRowCopy(int BlockRow, int MaxNumBlockEntries, 
00569            int & RowDim, int & NumBlockEntries, 
00570            int * BlockIndices, int * ColDims) const;
00571 
00573 
00588     int ExtractEntryCopy(int SizeOfValues, double * Values, int LDA, bool SumInto) const;
00589 
00591 
00603     int BeginExtractGlobalBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00604                int * & BlockIndices) const;
00605 
00607 
00619     int BeginExtractMyBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00620            int * & BlockIndices) const;
00621 
00622 
00624 
00631     int ExtractEntryView(Epetra_SerialDenseMatrix* & entry) const;
00632 
00634 
00648     int ExtractGlobalBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00649           int * & BlockIndices,
00650           Epetra_SerialDenseMatrix** & Values) const;
00651 
00653 
00667     int ExtractMyBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00668             int * & BlockIndices,
00669             Epetra_SerialDenseMatrix** & Values) const;
00670 
00671 
00673 
00679     int ExtractDiagonalCopy(Epetra_Vector & Diagonal) const;
00680 
00682 
00692     int BeginExtractBlockDiagonalCopy(int MaxNumBlockDiagonalEntries, 
00693               int & NumBlockDiagonalEntries, int * RowColDims ) const;
00695 
00710     int ExtractBlockDiagonalEntryCopy(int SizeOfValues, double * Values, int LDA, bool SumInto) const;
00711 
00713 
00721     int BeginExtractBlockDiagonalView(int & NumBlockDiagonalEntries, int * & RowColDims ) const;
00722 
00724 
00734     int ExtractBlockDiagonalEntryView(double * & Values, int & LDA) const;
00736 
00738 
00739 
00740 
00742 
00752     int Multiply1(bool TransA, const Epetra_Vector& x, Epetra_Vector& y) const;
00753 
00755 
00765     int Multiply(bool TransA, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00766 
00768 
00782     int Solve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_Vector& x, Epetra_Vector& y) const;
00783 
00785 
00799     int Solve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00800 
00801 
00803 
00813     int InvRowSums(Epetra_Vector& x) const;
00814 
00816 
00823     int LeftScale(const Epetra_Vector& x);
00824 
00826 
00836     int InvColSums(Epetra_Vector& x) const ;
00837 
00839 
00846     int RightScale(const Epetra_Vector& x);
00848 
00850 
00851 
00853     int OptimizeStorage();
00854 
00856     bool StorageOptimized() const {return(StorageOptimized_);};
00857 
00859     bool IndicesAreGlobal() const {return(Graph_->IndicesAreGlobal());};
00860 
00862     bool IndicesAreLocal() const {return(Graph_->IndicesAreLocal());};
00863 
00865     bool IndicesAreContiguous() const {return(Graph_->IndicesAreContiguous());};
00866 
00868     bool LowerTriangular() const {return(Graph_->LowerTriangular());};
00869 
00871     bool UpperTriangular() const {return(Graph_->UpperTriangular());};
00872 
00874     bool NoDiagonal() const {return(Graph_->NoDiagonal());};
00875 
00877   
00879 
00880 
00882     /* Returns the quantity \f$ \| A \|_\infty\f$ such that
00883        \f[\| A \|_\infty = \max_{1\lei\lem} \sum_{j=1}^n |a_{ij}| \f].
00884      \warning The NormInf() method will not properly calculate the infinity norm for a matrix that has entries that are
00885      replicated on multiple processors. */ 
00886     double NormInf() const;
00887 
00889     /* Returns the quantity \f$ \| A \|_1\f$ such that
00890        \f[\| A \|_1 = \max_{1\lej\len} \sum_{i=1}^m |a_{ij}| \f].
00891      \warning The NormOne() method will not properly calculate the one norm for a matrix that has entries that are
00892     */ 
00893     double NormOne() const;
00894 
00896     /* Returns the quantity \f[ \| A \|_{Frobenius} = \sqrt{\sum_{i=1}^m \sum_{j=1}^n\|a_{ij}\|^2}\f]
00897        \warning the NormFrobenius() method will not properly calculate the frobenius norm for a matrix that
00898        has entries which are replicated on multiple processors. In that case, the returned
00899        norm will be larger than the true norm.
00900     */
00901     double NormFrobenius() const;
00902 
00904     int MaxRowDim() const {return(Graph_->MaxRowDim());};
00905 
00907     int MaxColDim() const {return(Graph_->MaxColDim());};
00908 
00910     int GlobalMaxRowDim() const {return(Graph_->GlobalMaxRowDim());};
00911 
00913     int GlobalMaxColDim() const {return(Graph_->GlobalMaxColDim());};
00914 
00916     int NumMyRows() const {return(Graph_->NumMyRows());};
00918     int NumMyCols() const {return(Graph_->NumMyCols());};
00919 
00921     int NumMyNonzeros() const {return(Graph_->NumMyNonzeros());};
00922 
00924     int NumGlobalRows() const {return(Graph_->NumGlobalRows());};
00925 
00927     int NumGlobalCols() const {return(Graph_->NumGlobalCols());};
00928 
00930     /*
00931      Note that if maps are defined such that some nonzeros appear on
00932      multiple processors, then those nonzeros will be counted multiple
00933      times. If the user wishes to assemble a matrix from overlapping
00934      submatrices, they can use Epetra_FEVbrMatrix.
00935     */
00936     int NumGlobalNonzeros() const {return(Graph_->NumGlobalNonzeros());};
00937 
00939     int NumMyBlockRows() const {return(Graph_->NumMyBlockRows());};
00940 
00942     int NumMyBlockCols() const {return(Graph_->NumMyBlockCols());};
00943     
00945     int NumMyBlockEntries() const {return(Graph_->NumMyEntries());};
00946 
00948     int NumMyBlockDiagonals() const {return(Graph_->NumMyBlockDiagonals());};
00949     
00951     int NumMyDiagonals() const {return(Graph_->NumMyDiagonals());};
00952     
00954     int NumGlobalBlockRows() const {return(Graph_->NumGlobalBlockRows());};
00955     
00957     int NumGlobalBlockCols() const {return(Graph_->NumGlobalBlockCols());};
00958     
00960     int NumGlobalBlockEntries() const {return(Graph_->NumGlobalEntries());};
00961     
00963     int NumGlobalBlockDiagonals() const {return(Graph_->NumGlobalBlockDiagonals());};
00964     
00966     int NumGlobalDiagonals() const {return(Graph_->NumGlobalDiagonals());};
00967 
00969     int NumGlobalBlockEntries(int Row) const {return(Graph_->NumGlobalIndices(Row));};
00970 
00972     int NumAllocatedGlobalBlockEntries(int Row) const{return(Graph_->NumAllocatedGlobalIndices(Row));};
00973 
00975     int MaxNumBlockEntries() const {return(Graph_->MaxNumIndices());};
00976 
00978     int GlobalMaxNumBlockEntries() const {return(Graph_->GlobalMaxNumIndices());};
00979 
00981     int NumMyBlockEntries(int Row) const {return(Graph_->NumMyIndices(Row));};
00982 
00984     int NumAllocatedMyBlockEntries(int Row) const {return(Graph_->NumAllocatedMyIndices(Row));};
00985 
00987 
00991     int MaxNumNonzeros() const {return(Graph_->MaxNumNonzeros());};
00992 
00994 
00996     int GlobalMaxNumNonzeros() const {return(Graph_->GlobalMaxNumNonzeros());};
00997 
00999     int IndexBase() const {return(Graph_->IndexBase());};
01000 
01002     const Epetra_CrsGraph & Graph() const {return(*Graph_);};
01003 
01005     const Epetra_Import * Importer() const {return(Graph_->Importer());};
01006 
01008     const Epetra_Export * Exporter() const {return(Graph_->Exporter());};
01009 
01011     const Epetra_BlockMap & DomainMap() const {return(Graph_->DomainMap());};
01012 
01014     const Epetra_BlockMap & RangeMap() const  {return(Graph_->RangeMap());};
01015 
01017     const Epetra_BlockMap & RowMap() const {return(Graph_->RowMap());};
01018 
01020     const Epetra_BlockMap & ColMap() const {return(Graph_->ColMap());};
01021 
01023 
01025     const Epetra_Comm & Comm() const {return(Graph_->Comm());};
01026 
01028   
01030 
01031 
01032     int LRID( int GRID_in) const {return(Graph_->LRID(GRID_in));};
01033 
01035     int GRID( int LRID_in) const {return(Graph_->GRID(LRID_in));};
01036 
01038     int LCID( int GCID_in) const {return(Graph_->LCID(GCID_in));};
01039 
01041     int GCID( int LCID_in) const {return(Graph_->GCID(LCID_in));};
01042  
01044     bool  MyGRID(int GRID_in) const {return(Graph_->MyGRID(GRID_in));};
01045    
01047     bool  MyLRID(int LRID_in) const {return(Graph_->MyLRID(LRID_in));};
01048 
01050     bool  MyGCID(int GCID_in) const {return(Graph_->MyGCID(GCID_in));};
01051    
01053     bool  MyLCID(int LCID_in) const {return(Graph_->MyLCID(LCID_in));};
01054 
01056     bool MyGlobalBlockRow(int GID) const {return(Graph_->MyGlobalRow(GID));};
01058   
01060 
01061 
01063   virtual void Print(ostream & os) const;
01065 
01067 
01068 
01070     const char * Label() const {return(Epetra_Object::Label());};
01071     
01073 
01082   int SetUseTranspose(bool UseTranspose_in) {UseTranspose_ = UseTranspose_in; return(0);};
01083 
01085 
01093   int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
01094 
01096 
01109   int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
01110 
01112     bool HasNormInf() const {return(true);};
01113 
01115     bool UseTranspose() const {return(UseTranspose_);};
01116 
01118     const Epetra_Map & OperatorDomainMap() const 
01119     {
01120       if (!HavePointObjects_) GeneratePointObjects();
01121       if (UseTranspose()) return(*OperatorRangeMap_);
01122       else return(*OperatorDomainMap_);
01123     }
01124 
01126     const Epetra_Map & OperatorRangeMap() const 
01127     {
01128       if (!HavePointObjects_) GeneratePointObjects();
01129       if (UseTranspose()) return(*OperatorDomainMap_);
01130       else return(*OperatorRangeMap_);
01131     }
01132 
01134 
01135 
01136 
01138 
01152     int ExtractGlobalRowCopy(int GlobalRow, int Length, int & NumEntries, double *Values, int * Indices) const;
01153 
01155 
01169     int ExtractMyRowCopy(int MyRow, int Length, int & NumEntries, double *Values, int * Indices) const;
01170 
01172 
01180     int NumMyRowEntries(int MyRow, int & NumEntries) const;
01181 
01183     int MaxNumEntries() const;
01184 
01186     const Epetra_BlockMap& Map() const { return Epetra_DistObject::Map(); }
01187 
01189     const Epetra_Map & RowMatrixRowMap() const 
01190       { if (!HavePointObjects_) GeneratePointObjects(); return(*RowMatrixRowMap_); };
01191 
01193     const Epetra_Map & RowMatrixColMap() const 
01194       { if (!HavePointObjects_) GeneratePointObjects(); return(*RowMatrixColMap_); };
01195 
01197     const Epetra_Import * RowMatrixImporter() const 
01198       { if (!HavePointObjects_) GeneratePointObjects(); return(RowMatrixImporter_); };
01199 
01201 
01203 
01204 
01206     const Epetra_BlockMap & BlockImportMap() const {return(Graph_->ImportMap());};
01207 
01209     int TransformToLocal();
01210     
01212     int TransformToLocal(const Epetra_BlockMap* DomainMap, const Epetra_BlockMap* RangeMap);
01214 
01215  protected:
01216     void DeleteMemory();
01217     bool Allocated() const {return(Allocated_);};
01218     int SetAllocated(bool Flag) {Allocated_ = Flag; return(0);};
01219     Epetra_SerialDenseMatrix *** Values() const {return(Entries_);};
01220 
01221   // Internal utilities
01222 
01223   int DoMultiply(bool TransA, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
01224   int DoSolve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
01225   void InitializeDefaults();
01226   int Allocate();
01227   int BeginInsertValues(int BlockRow, int NumBlockEntries, 
01228       int * BlockIndices, bool IndicesAreLocal);
01229   int BeginReplaceValues(int BlockRow, int NumBlockEntries, 
01230        int *BlockIndices, bool IndicesAreLocal);
01231   int BeginSumIntoValues(int BlockRow, int NumBlockEntries, 
01232        int *BlockIndices, bool IndicesAreLocal);
01233   int SetupForSubmits(int BlockRow, int NumBlockEntries, int * BlockIndices, 
01234           bool IndicesAreLocal, Epetra_CombineMode SubmitMode);
01235   int EndReplaceSumIntoValues();
01236   int EndInsertValues();
01237 
01238   int CopyMat(double * A, int LDA, int NumRows, int NumCols, 
01239         double * B, int LDB, bool SumInto) const;
01240   int BeginExtractBlockRowCopy(int BlockRow, int MaxNumBlockEntries, 
01241              int & RowDim, int & NumBlockEntries, 
01242              int * BlockIndices, int * ColDims, 
01243              bool IndicesAreLocal) const;
01244   int SetupForExtracts(int BlockRow, int & RowDim, int NumBlockEntries,
01245            bool ExtractView, bool IndicesAreLocal) const;
01246   int ExtractBlockDimsCopy(int NumBlockEntries, int * ColDims) const;
01247   int ExtractBlockRowPointers(int BlockRow, int MaxNumBlockEntries, 
01248             int & RowDim, int & NumBlockEntries, 
01249             int * BlockIndices, 
01250             Epetra_SerialDenseMatrix ** & Values,
01251             bool IndicesAreLocal) const;
01252   int BeginExtractBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
01253              int * & BlockIndices, 
01254              bool IndicesAreLocal) const;
01255   int CopyMatDiag(double * A, int LDA, int NumRows, int NumCols, 
01256       double * Diagonal) const;
01257   int ReplaceMatDiag(double * A, int LDA, int NumRows, int NumCols, 
01258          double * Diagonal);
01259 
01260   //This BlockRowMultiply accepts Alpha and Beta arguments. It is called
01261   //from within the 'solve' methods.
01262   void BlockRowMultiply(bool TransA, int RowDim, int NumEntries, 
01263       int * BlockIndices, int RowOff,
01264       int * FirstPointInElementList, int * ElementSizeList,
01265       double Alpha, Epetra_SerialDenseMatrix** As,
01266       double ** X, double Beta, double ** Y, int NumVectors) const;
01267 
01268   //This BlockRowMultiply doesn't accept Alpha and Beta arguments, instead it
01269   //assumes that they are both 1.0. It is called from within the 'Multiply'
01270   //methods.
01271   void BlockRowMultiply(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   // Assumes Alpha=Beta=1 and works only on storage optimized matrices
01279   //
01280   void FastBlockRowMultiply(bool TransA, int RowDim, int NumEntries, 
01281           int * BlockIndices, int RowOff,
01282           int * FirstPointInElementList,
01283           int * ElementSizeList,
01284           Epetra_SerialDenseMatrix** As,
01285           double ** X, double ** Y, int NumVectors) const;
01286 
01287   int InverseSums(bool DoRows, Epetra_Vector& x) const;
01288   int Scale(bool DoRows, const Epetra_Vector& x);
01289   void BlockRowNormInf(int RowDim, int NumEntries, 
01290            Epetra_SerialDenseMatrix** As, 
01291            double * Y) const;
01292   void BlockRowNormOne(int RowDim, int NumEntries, int * BlockRowIndices,
01293            Epetra_SerialDenseMatrix** As, 
01294            int * ColFirstPointInElementList, double * x) const;
01295   void SetStaticGraph(bool Flag) {StaticGraph_ = Flag;};
01296 
01297   int CheckSizes(const Epetra_SrcDistObject& A);
01298 
01299   int CopyAndPermute(const Epetra_SrcDistObject & Source,
01300                      int NumSameIDs, 
01301                      int NumPermuteIDs,
01302                      int * PermuteToLIDs,
01303                      int *PermuteFromLIDs,
01304                      const Epetra_OffsetIndex * Indexor);
01305   
01306   int PackAndPrepare(const Epetra_SrcDistObject & Source,
01307                      int NumExportIDs,
01308                      int * ExportLIDs,
01309                      int & LenExports,
01310                      char * & Exports,
01311                      int & SizeOfPacket,
01312                      int * Sizes,
01313                      bool & VarSizes,
01314                      Epetra_Distributor & Distor);
01315   
01316   int UnpackAndCombine(const Epetra_SrcDistObject & Source, 
01317                        int NumImportIDs,
01318                        int * ImportLIDs, 
01319                        int LenImports,
01320                        char * Imports,
01321                        int & SizeOfPacket, 
01322                        Epetra_Distributor & Distor,
01323                        Epetra_CombineMode CombineMode,
01324                        const Epetra_OffsetIndex * Indexor);
01325 
01327   int SortEntries();
01328 
01330   bool Sorted() const {return(Graph_->Sorted());};
01331 
01333   int MergeRedundantEntries();
01334 
01336   bool NoRedundancies() const {return(Graph_->NoRedundancies());};
01337 
01338   bool StaticGraph() const {return(StaticGraph_);};
01339 
01340   int GeneratePointObjects() const;
01341   int BlockMap2PointMap(const Epetra_BlockMap & BlockMap, Epetra_Map * & PointMap) const;
01342   int UpdateOperatorXY(const Epetra_MultiVector& X, const Epetra_MultiVector& Y) const;
01343 
01344   Epetra_CrsGraph * Graph_;
01345   bool Allocated_;
01346   bool StaticGraph_;
01347   bool UseTranspose_;
01348   bool constructedWithFilledGraph_;
01349   bool matrixFillCompleteCalled_;
01350   bool StorageOptimized_;
01351 
01352   int NumMyBlockRows_;
01353 
01354   Epetra_DataAccess CV_;
01355 
01356 
01357   int * NumBlockEntriesPerRow_;
01358   int * NumAllocatedBlockEntriesPerRow_;
01359   int ** Indices_;
01360   int * ElementSizeList_;
01361   int * FirstPointInElementList_;
01362 
01363   Epetra_SerialDenseMatrix ***Entries_;
01364 
01365   double *All_Values_Orig_;
01366   double *All_Values_;
01367 
01368   mutable double NormInf_;
01369   mutable double NormOne_;
01370   mutable double NormFrob_;
01371 
01372   mutable Epetra_MultiVector * ImportVector_;
01373   mutable Epetra_MultiVector * ExportVector_;
01374 
01375   // State variables needed for constructing matrix entry-by-entry
01376   mutable int *TempRowDims_;
01377   mutable Epetra_SerialDenseMatrix **TempEntries_;
01378   mutable int LenTemps_;
01379   mutable int CurBlockRow_;
01380   mutable int CurNumBlockEntries_;
01381   mutable int * CurBlockIndices_;
01382   mutable int CurEntry_;
01383   mutable bool CurIndicesAreLocal_;
01384   mutable Epetra_CombineMode CurSubmitMode_;
01385   
01386   // State variables needed for extracting entries
01387   mutable int CurExtractBlockRow_;
01388   mutable int CurExtractEntry_; 
01389   mutable int CurExtractNumBlockEntries_;
01390   mutable bool CurExtractIndicesAreLocal_;
01391   mutable bool CurExtractView_;
01392   mutable int CurRowDim_;
01393 
01394   // State variable for extracting block diagonal entries
01395   mutable int CurBlockDiag_;
01396 
01397   // Maps and importer that support the Epetra_RowMatrix interface
01398   mutable Epetra_Map * RowMatrixRowMap_;
01399   mutable Epetra_Map * RowMatrixColMap_;
01400   mutable Epetra_Import * RowMatrixImporter_;
01401 
01402   // Maps that support the Epetra_Operator interface
01403   mutable Epetra_Map * OperatorDomainMap_;
01404   mutable Epetra_Map * OperatorRangeMap_;
01405   mutable Epetra_MultiVector * OperatorX_;
01406   mutable Epetra_MultiVector * OperatorY_;
01407 
01408   // bool to indicate if above four point maps and importer have already been created
01409   mutable bool HavePointObjects_;
01410 
01411   bool squareFillCompleteCalled_;
01412 };
01413 
01414 #endif /* EPETRA_VBRMATRIX_H */
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines