Epetra_VbrMatrix.h

Go to the documentation of this file.
00001 //@HEADER
00002 /*
00003 ************************************************************************
00004 
00005               Epetra: Linear Algebra Services Package 
00006                 Copyright (2001) Sandia Corporation
00007 
00008 Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00009 license for use of this work by or on behalf of the U.S. Government.
00010 
00011 This library is free software; you can redistribute it and/or modify
00012 it under the terms of the GNU Lesser General Public License as
00013 published by the Free Software Foundation; either version 2.1 of the
00014 License, or (at your option) any later version.
00015  
00016 This library is distributed in the hope that it will be useful, but
00017 WITHOUT ANY WARRANTY; without even the implied warranty of
00018 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00019 Lesser General Public License for more details.
00020  
00021 You should have received a copy of the GNU Lesser General Public
00022 License along with this library; if not, write to the Free Software
00023 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00024 USA
00025 Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00026 
00027 ************************************************************************
00028 */
00029 //@HEADER
00030 
00031 #ifndef EPETRA_VBRMATRIX_H
00032 #define EPETRA_VBRMATRIX_H
00033 
00034 #include <Epetra_DistObject.h> 
00035 #include <Epetra_CompObject.h> 
00036 #include <Epetra_BLAS.h>
00037 #include <Epetra_RowMatrix.h>
00038 #include <Epetra_Operator.h>
00039 #include <Epetra_CrsGraph.h>
00040 #include <Epetra_SerialDenseMatrix.h>
00041 class Epetra_BlockMap;
00042 class Epetra_Map;
00043 class Epetra_Import;
00044 class Epetra_Export;
00045 class Epetra_Vector;
00046 class Epetra_MultiVector;
00047 
00049 
00147 class Epetra_VbrMatrix : public Epetra_DistObject,
00148           public Epetra_CompObject,
00149           public Epetra_BLAS,
00150           public virtual Epetra_RowMatrix {
00151  public:
00152 
00154 
00155 
00156 
00166   Epetra_VbrMatrix(Epetra_DataAccess CV, const Epetra_BlockMap& RowMap, int *NumBlockEntriesPerRow);
00167   
00169 
00180   Epetra_VbrMatrix(Epetra_DataAccess CV, const Epetra_BlockMap& RowMap, int NumBlockEntriesPerRow);
00181 
00183 
00195   Epetra_VbrMatrix(Epetra_DataAccess CV, const Epetra_BlockMap& RowMap, const Epetra_BlockMap& ColMap, int *NumBlockEntriesPerRow);
00196   
00198 
00211   Epetra_VbrMatrix(Epetra_DataAccess CV, const Epetra_BlockMap& RowMap, const Epetra_BlockMap& ColMap, int NumBlockEntriesPerRow);
00212 
00214 
00223   Epetra_VbrMatrix(Epetra_DataAccess CV, const Epetra_CrsGraph & Graph);
00224 
00226   Epetra_VbrMatrix(const Epetra_VbrMatrix & Matrix);
00227 
00229   virtual ~Epetra_VbrMatrix();
00231   
00233 
00234 
00235   Epetra_VbrMatrix& operator=(const Epetra_VbrMatrix& src);
00236 
00238 
00244     int PutScalar(double ScalarConstant);
00245 
00247 
00253     int Scale(double ScalarConstant);
00254 
00255 
00257 
00267     int BeginInsertGlobalValues(int BlockRow,
00268         int NumBlockEntries,
00269         int * BlockIndices);
00270 
00272 
00282     int BeginInsertMyValues(int BlockRow, int NumBlockEntries, int * BlockIndices);
00283 
00285 
00295     int BeginReplaceGlobalValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00296 
00298 
00308     int BeginReplaceMyValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00309 
00311 
00321     int BeginSumIntoGlobalValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00322 
00324 
00334     int BeginSumIntoMyValues(int BlockRow, int NumBlockEntries, int *BlockIndices);
00335 
00337     /* Submit a block entry that will recorded in the block row that was initiated by one of the
00338        Begin routines listed above.  Once a one of the following routines: BeginInsertGlobalValues(),
00339        BeginInsertMyValues(), BeginReplaceGlobalValues(), BeginReplaceMyValues(), BeginSumIntoGlobalValues(),
00340        BeginSumIntoMyValues(), you \e must call SubmitBlockEntry() NumBlockEntries times to register the values 
00341        corresponding to the block indices passed in to the Begin routine.  If the Epetra_VbrMatrix constuctor
00342        was called in Copy mode, the values will be copied.  However, no copying will be done until the EndSubmitEntries()
00343        function is call to complete submission of the current block row.  If the constructor was called in View mode, all
00344        block entries passed via SubmitBlockEntry() will not be copied, but a pointer will be set to point to the argument Values
00345        that was passed in by the user.
00346 
00347        For performance reasons, SubmitBlockEntry() does minimal processing of data.  Any processing that can be
00348        delayed is performed in EndSubmitEntries().
00349 
00350     \param In
00351            Values - The starting address of the values.
00352     \param In
00353            LDA - The stride between successive columns of Values.
00354     \param In
00355            NumRows - The number of rows passed in.
00356     \param In
00357            NumCols - The number of columns passed in.
00358 
00359     \return Integer error code, set to 0 if successful.
00360     */
00361     int SubmitBlockEntry(double *Values, int LDA, int NumRows, int NumCols);
00362 
00364     /* Submit a block entry that will recorded in the block row that was initiated by one of the
00365        Begin routines listed above.  Once a one of the following routines: BeginInsertGlobalValues(),
00366        BeginInsertMyValues(), BeginReplaceGlobalValues(), BeginReplaceMyValues(), BeginSumIntoGlobalValues(),
00367        BeginSumIntoMyValues(), you \e must call SubmitBlockEntry() NumBlockEntries times to register the values 
00368        corresponding to the block indices passed in to the Begin routine.  If the Epetra_VbrMatrix constuctor
00369        was called in Copy mode, the values will be copied.  However, no copying will be done until the EndSubmitEntries()
00370        function is call to complete submission of the current block row.  If the constructor was called in View mode, all
00371        block entries passed via SubmitBlockEntry() will not be copied, but a pointer will be set to point to the argument Values
00372        that was passed in by the user.
00373 
00374        For performance reasons, SubmitBlockEntry() does minimal processing of data.  Any processing that can be
00375        delayed is performed in EndSubmitEntries().
00376 
00377     \param In
00378            Mat - Preformed dense matrix block.
00379 
00380     \return Integer error code, set to 0 if successful.
00381     */
00382     int SubmitBlockEntry( Epetra_SerialDenseMatrix &Mat );
00383 
00385 
00390     int EndSubmitEntries();
00391 
00393 
00404     int ReplaceDiagonalValues(const Epetra_Vector & Diagonal);
00405 
00407     /* This version of FillComplete assumes that the domain and range
00408        distributions are identical to the matrix row distributions.
00409        \return error code, 0 if successful. Returns a positive warning code of 3
00410         if the matrix is rectangular (meaning that the other overloading of
00411         FillComplete should have been called, with differen domain-map and
00412         range-map specified).
00413     */
00414     int FillComplete();
00415 
00417     /* This version of FillComplete requires the explicit specification of the domain
00418        and range distribution maps.  These maps are used for importing and exporting vector
00419        and multi-vector elements that are needed for distributed matrix computations.  For
00420        example, to compute y = Ax in parallel, we would specify the DomainMap as the distribution
00421        of the vector x and the RangeMap as the distribution of the vector y.
00422     \param In
00423            DomainMap - Map that describes the distribution of vector and multi-vectors in the
00424                  matrix domain.
00425     \param In
00426            RangeMap - Map that describes the distribution of vector and multi-vectors in the
00427                  matrix range.
00428 
00429     \return error code, 0 if successful. positive warning code of 2 if it is detected that the
00430     matrix-graph got out of sync since this matrix was constructed (for instance if
00431     graph.FillComplete() was called by another matrix that shares the graph)
00432     */
00433     int FillComplete(const Epetra_BlockMap& DomainMap, const Epetra_BlockMap& RangeMap);
00434 
00436     bool Filled() const {return(Graph_->Filled());};
00438 
00440 
00441 
00443 
00465     int ExtractGlobalBlockRowPointers(int BlockRow, int MaxNumBlockEntries, 
00466               int & RowDim,  int & NumBlockEntries, 
00467               int * BlockIndices,
00468               Epetra_SerialDenseMatrix ** & Values) const;
00469 
00471 
00493     int ExtractMyBlockRowPointers(int BlockRow, int MaxNumBlockEntries, 
00494           int & RowDim, int & NumBlockEntries, 
00495           int * BlockIndices,
00496           Epetra_SerialDenseMatrix** & Values) const;
00497 
00499 
00515     int BeginExtractGlobalBlockRowCopy(int BlockRow, int MaxNumBlockEntries, 
00516                int & RowDim,  int & NumBlockEntries, 
00517                int * BlockIndices, int * ColDims) const;
00518 
00520 
00536     int BeginExtractMyBlockRowCopy(int BlockRow, int MaxNumBlockEntries, 
00537            int & RowDim, int & NumBlockEntries, 
00538            int * BlockIndices, int * ColDims) const;
00539 
00541 
00556     int ExtractEntryCopy(int SizeOfValues, double * Values, int LDA, bool SumInto) const;
00557 
00559 
00571     int BeginExtractGlobalBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00572                int * & BlockIndices) const;
00573 
00575 
00587     int BeginExtractMyBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00588            int * & BlockIndices) const;
00589 
00590 
00592 
00599     int ExtractEntryView(Epetra_SerialDenseMatrix* & entry) const;
00600 
00602 
00616     int ExtractGlobalBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00617           int * & BlockIndices,
00618           Epetra_SerialDenseMatrix** & Values) const;
00619 
00621 
00635     int ExtractMyBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
00636             int * & BlockIndices,
00637             Epetra_SerialDenseMatrix** & Values) const;
00638 
00639 
00641 
00647     int ExtractDiagonalCopy(Epetra_Vector & Diagonal) const;
00648 
00650 
00660     int BeginExtractBlockDiagonalCopy(int MaxNumBlockDiagonalEntries, 
00661               int & NumBlockDiagonalEntries, int * RowColDims ) const;
00663 
00678     int ExtractBlockDiagonalEntryCopy(int SizeOfValues, double * Values, int LDA, bool SumInto) const;
00679 
00681 
00689     int BeginExtractBlockDiagonalView(int & NumBlockDiagonalEntries, int * & RowColDims ) const;
00690 
00692 
00702     int ExtractBlockDiagonalEntryView(double * & Values, int & LDA) const;
00704 
00706 
00707 
00708 
00710 
00720     int Multiply1(bool TransA, const Epetra_Vector& x, Epetra_Vector& y) const;
00721 
00723 
00733     int Multiply(bool TransA, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00734 
00736 
00750     int Solve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_Vector& x, Epetra_Vector& y) const;
00751 
00753 
00767     int Solve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00768 
00769 
00771 
00781     int InvRowSums(Epetra_Vector& x) const;
00782 
00784 
00791     int LeftScale(const Epetra_Vector& x);
00792 
00794 
00804     int InvColSums(Epetra_Vector& x) const ;
00805 
00807 
00814     int RightScale(const Epetra_Vector& x);
00816 
00818 
00819 
00821     int OptimizeStorage();
00822 
00824     bool StorageOptimized() const {return(StorageOptimized_);};
00825 
00827     bool IndicesAreGlobal() const {return(Graph_->IndicesAreGlobal());};
00828 
00830     bool IndicesAreLocal() const {return(Graph_->IndicesAreLocal());};
00831 
00833     bool IndicesAreContiguous() const {return(Graph_->IndicesAreContiguous());};
00834 
00836     bool LowerTriangular() const {return(Graph_->LowerTriangular());};
00837 
00839     bool UpperTriangular() const {return(Graph_->UpperTriangular());};
00840 
00842     bool NoDiagonal() const {return(Graph_->NoDiagonal());};
00843 
00845   
00847 
00848 
00850     /* Returns the quantity \f$ \| A \|_\infty\f$ such that
00851        \f[\| A \|_\infty = \max_{1\lei\lem} \sum_{j=1}^n |a_{ij}| \f].
00852      \warning The NormInf() method will not properly calculate the infinity norm for a matrix that has entries that are
00853      replicated on multiple processors. */ 
00854     double NormInf() const;
00855 
00857     /* Returns the quantity \f$ \| A \|_1\f$ such that
00858        \f[\| A \|_1 = \max_{1\lej\len} \sum_{i=1}^m |a_{ij}| \f].
00859      \warning The NormOne() method will not properly calculate the one norm for a matrix that has entries that are
00860     */ 
00861     double NormOne() const;
00862 
00864     /* Returns the quantity \f[ \| A \|_{Frobenius} = \sqrt{\sum_{i=1}^m \sum_{j=1}^n\|a_{ij}\|^2}\f]
00865        \warning the NormFrobenius() method will not properly calculate the frobenius norm for a matrix that
00866        has entries which are replicated on multiple processors. In that case, the returned
00867        norm will be larger than the true norm.
00868     */
00869     double NormFrobenius() const;
00870 
00872     int MaxRowDim() const {return(Graph_->MaxRowDim());};
00873 
00875     int MaxColDim() const {return(Graph_->MaxColDim());};
00876 
00878     int GlobalMaxRowDim() const {return(Graph_->GlobalMaxRowDim());};
00879 
00881     int GlobalMaxColDim() const {return(Graph_->GlobalMaxColDim());};
00882 
00884     int NumMyRows() const {return(Graph_->NumMyRows());};
00886     int NumMyCols() const {return(Graph_->NumMyCols());};
00887 
00889     int NumMyNonzeros() const {return(Graph_->NumMyNonzeros());};
00890 
00892     int NumGlobalRows() const {return(Graph_->NumGlobalRows());};
00893 
00895     int NumGlobalCols() const {return(Graph_->NumGlobalCols());};
00896 
00898     /*
00899      Note that if maps are defined such that some nonzeros appear on
00900      multiple processors, then those nonzeros will be counted multiple
00901      times. If the user wishes to assemble a matrix from overlapping
00902      submatrices, they can use Epetra_FEVbrMatrix.
00903     */
00904     int NumGlobalNonzeros() const {return(Graph_->NumGlobalNonzeros());};
00905 
00907     int NumMyBlockRows() const {return(Graph_->NumMyBlockRows());};
00908 
00910     int NumMyBlockCols() const {return(Graph_->NumMyBlockCols());};
00911     
00913     int NumMyBlockEntries() const {return(Graph_->NumMyEntries());};
00914 
00916     int NumMyBlockDiagonals() const {return(Graph_->NumMyBlockDiagonals());};
00917     
00919     int NumMyDiagonals() const {return(Graph_->NumMyDiagonals());};
00920     
00922     int NumGlobalBlockRows() const {return(Graph_->NumGlobalBlockRows());};
00923     
00925     int NumGlobalBlockCols() const {return(Graph_->NumGlobalBlockCols());};
00926     
00928     int NumGlobalBlockEntries() const {return(Graph_->NumGlobalEntries());};
00929     
00931     int NumGlobalBlockDiagonals() const {return(Graph_->NumGlobalBlockDiagonals());};
00932     
00934     int NumGlobalDiagonals() const {return(Graph_->NumGlobalDiagonals());};
00935 
00937     int NumGlobalBlockEntries(int Row) const {return(Graph_->NumGlobalIndices(Row));};
00938 
00940     int NumAllocatedGlobalBlockEntries(int Row) const{return(Graph_->NumAllocatedGlobalIndices(Row));};
00941 
00943     int MaxNumBlockEntries() const {return(Graph_->MaxNumIndices());};
00944 
00946     int GlobalMaxNumBlockEntries() const {return(Graph_->GlobalMaxNumIndices());};
00947 
00949     int NumMyBlockEntries(int Row) const {return(Graph_->NumMyIndices(Row));};
00950 
00952     int NumAllocatedMyBlockEntries(int Row) const {return(Graph_->NumAllocatedMyIndices(Row));};
00953 
00955 
00959     int MaxNumNonzeros() const {return(Graph_->MaxNumNonzeros());};
00960 
00962 
00964     int GlobalMaxNumNonzeros() const {return(Graph_->GlobalMaxNumNonzeros());};
00965 
00967     int IndexBase() const {return(Graph_->IndexBase());};
00968 
00970     const Epetra_CrsGraph & Graph() const {return(*Graph_);};
00971 
00973     const Epetra_Import * Importer() const {return(Graph_->Importer());};
00974 
00976     const Epetra_Export * Exporter() const {return(Graph_->Exporter());};
00977 
00979     const Epetra_BlockMap & DomainMap() const {return(Graph_->DomainMap());};
00980 
00982     const Epetra_BlockMap & RangeMap() const  {return(Graph_->RangeMap());};
00983 
00985     const Epetra_BlockMap & RowMap() const {return(Graph_->RowMap());};
00986 
00988     const Epetra_BlockMap & ColMap() const {return(Graph_->ColMap());};
00989 
00991 
00993     const Epetra_Comm & Comm() const {return(Graph_->Comm());};
00994 
00996   
00998 
00999 
01000     int LRID( int GRID) const {return(Graph_->LRID(GRID));};
01001 
01003     int GRID( int LRID) const {return(Graph_->GRID(LRID));};
01004 
01006     int LCID( int GCID) const {return(Graph_->LCID(GCID));};
01007 
01009     int GCID( int LCID) const {return(Graph_->GCID(LCID));};
01010  
01012     bool  MyGRID(int GRID) const {return(Graph_->MyGRID(GRID));};
01013    
01015     bool  MyLRID(int LRID) const {return(Graph_->MyLRID(LRID));};
01016 
01018     bool  MyGCID(int GCID) const {return(Graph_->MyGCID(GCID));};
01019    
01021     bool  MyLCID(int LCID) const {return(Graph_->MyLCID(LCID));};
01022 
01024     bool MyGlobalBlockRow(int GID) const {return(Graph_->MyGlobalRow(GID));};
01026   
01028 
01029 
01031   virtual void Print(ostream & os) const;
01033 
01035 
01036 
01038     const char * Label() const {return(Epetra_Object::Label());};
01039     
01041 
01050   int SetUseTranspose(bool UseTranspose) {UseTranspose_ = UseTranspose; return(0);};
01051 
01053 
01061   int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
01062 
01064 
01077   int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
01078 
01080     bool HasNormInf() const {return(true);};
01081 
01083     bool UseTranspose() const {return(UseTranspose_);};
01084 
01086     const Epetra_Map & OperatorDomainMap() const 
01087     {
01088       if (!HavePointObjects_) GeneratePointObjects();
01089       if (UseTranspose()) return(*OperatorRangeMap_);
01090       else return(*OperatorDomainMap_);
01091     }
01092 
01094     const Epetra_Map & OperatorRangeMap() const 
01095     {
01096       if (!HavePointObjects_) GeneratePointObjects();
01097       if (UseTranspose()) return(*OperatorDomainMap_);
01098       else return(*OperatorRangeMap_);
01099     }
01100 
01102 
01103 
01104 
01106 
01120     int ExtractGlobalRowCopy(int GlobalRow, int Length, int & NumEntries, double *Values, int * Indices) const;
01121 
01123 
01137     int ExtractMyRowCopy(int MyRow, int Length, int & NumEntries, double *Values, int * Indices) const;
01138 
01140 
01148     int NumMyRowEntries(int MyRow, int & NumEntries) const;
01149 
01151     int MaxNumEntries() const;
01152 
01154     const Epetra_Map & RowMatrixRowMap() const 
01155       { if (!HavePointObjects_) GeneratePointObjects(); return(*RowMatrixRowMap_); };
01156 
01158     const Epetra_Map & RowMatrixColMap() const 
01159       { if (!HavePointObjects_) GeneratePointObjects(); return(*RowMatrixColMap_); };
01160 
01162     const Epetra_Import * RowMatrixImporter() const 
01163       { if (!HavePointObjects_) GeneratePointObjects(); return(RowMatrixImporter_); };
01164 
01166 
01168 
01169 
01171     const Epetra_BlockMap & BlockImportMap() const {return(Graph_->ImportMap());};
01172 
01174     int TransformToLocal();
01175     
01177     int TransformToLocal(const Epetra_BlockMap* DomainMap, const Epetra_BlockMap* RangeMap);
01179 
01180  protected:
01181     void DeleteMemory();
01182     bool Allocated() const {return(Allocated_);};
01183     int SetAllocated(bool Flag) {Allocated_ = Flag; return(0);};
01184     Epetra_SerialDenseMatrix *** Values() const {return(Entries_);};
01185 
01186   // Internal utilities
01187 
01188   void InitializeDefaults();
01189   int Allocate();
01190   int BeginInsertValues(int BlockRow, int NumBlockEntries, 
01191       int * BlockIndices, bool IndicesAreLocal);
01192   int BeginReplaceValues(int BlockRow, int NumBlockEntries, 
01193        int *BlockIndices, bool IndicesAreLocal);
01194   int BeginSumIntoValues(int BlockRow, int NumBlockEntries, 
01195        int *BlockIndices, bool IndicesAreLocal);
01196   int SetupForSubmits(int BlockRow, int NumBlockEntries, int * BlockIndices, 
01197           bool IndicesAreLocal, Epetra_CombineMode SubmitMode);
01198   int EndReplaceSumIntoValues();
01199   int EndInsertValues();
01200 
01201   int CopyMat(double * A, int LDA, int NumRows, int NumCols, 
01202         double * B, int LDB, bool SumInto) const;
01203   int BeginExtractBlockRowCopy(int BlockRow, int MaxNumBlockEntries, 
01204              int & RowDim, int & NumBlockEntries, 
01205              int * BlockIndices, int * ColDims, 
01206              bool IndicesAreLocal) const;
01207   int SetupForExtracts(int BlockRow, int & RowDim, int NumBlockEntries,
01208            bool ExtractView, bool IndicesAreLocal) const;
01209   int ExtractBlockDimsCopy(int NumBlockEntries, int * ColDims) const;
01210   int ExtractBlockRowPointers(int BlockRow, int MaxNumBlockEntries, 
01211             int & RowDim, int & NumBlockEntries, 
01212             int * BlockIndices, 
01213             Epetra_SerialDenseMatrix ** & Values,
01214             bool IndicesAreLocal) const;
01215   int BeginExtractBlockRowView(int BlockRow, int & RowDim, int & NumBlockEntries, 
01216              int * & BlockIndices, 
01217              bool IndicesAreLocal) const;
01218   int CopyMatDiag(double * A, int LDA, int NumRows, int NumCols, 
01219       double * Diagonal) const;
01220   int ReplaceMatDiag(double * A, int LDA, int NumRows, int NumCols, 
01221          double * Diagonal);
01222 
01223   //This BlockRowMultiply accepts Alpha and Beta arguments. It is called
01224   //from within the 'solve' methods.
01225   void BlockRowMultiply(bool TransA, int RowDim, int NumEntries, 
01226       int * BlockIndices, int RowOff,
01227       int * FirstPointInElementList, int * ElementSizeList,
01228       double Alpha, Epetra_SerialDenseMatrix** As,
01229       double ** X, double Beta, double ** Y, int NumVectors) const;
01230 
01231   //This BlockRowMultiply doesn't accept Alpha and Beta arguments, instead it
01232   //assumes that they are both 1.0. It is called from within the 'Multiply'
01233   //methods.
01234   void BlockRowMultiply(bool TransA, int RowDim, int NumEntries, 
01235       int * BlockIndices, int RowOff,
01236       int * FirstPointInElementList,
01237       int * ElementSizeList,
01238       Epetra_SerialDenseMatrix** As,
01239       double ** X, double ** Y, int NumVectors) const;
01240   //
01241   // Assumes Alpha=Beta=1 and works only on storage optimized matrices
01242   //
01243   void FastBlockRowMultiply(bool TransA, int RowDim, int NumEntries, 
01244           int * BlockIndices, int RowOff,
01245           int * FirstPointInElementList,
01246           int * ElementSizeList,
01247           Epetra_SerialDenseMatrix** As,
01248           double ** X, double ** Y, int NumVectors) const;
01249 
01250   int InverseSums(bool DoRows, Epetra_Vector& x) const;
01251   int Scale(bool DoRows, const Epetra_Vector& x);
01252   void BlockRowNormInf(int RowDim, int NumEntries, 
01253            Epetra_SerialDenseMatrix** As, 
01254            double * Y) const;
01255   void BlockRowNormOne(int RowDim, int NumEntries, int * BlockRowIndices,
01256            Epetra_SerialDenseMatrix** As, 
01257            int * ColFirstPointInElementList, double * x) const;
01258   void SetStaticGraph(bool Flag) {StaticGraph_ = Flag;};
01259 
01260   int CheckSizes(const Epetra_SrcDistObject& A);
01261 
01262   int CopyAndPermute(const Epetra_SrcDistObject & Source,
01263                      int NumSameIDs, 
01264                      int NumPermuteIDs,
01265                      int * PermuteToLIDs,
01266                      int *PermuteFromLIDs,
01267                      const Epetra_OffsetIndex * Indexor);
01268   
01269   int PackAndPrepare(const Epetra_SrcDistObject & Source,
01270                      int NumExportIDs,
01271                      int * ExportLIDs,
01272                      int & LenExports,
01273                      char * & Exports,
01274                      int & SizeOfPacket,
01275                      int * Sizes,
01276                      bool & VarSizes,
01277                      Epetra_Distributor & Distor);
01278   
01279   int UnpackAndCombine(const Epetra_SrcDistObject & Source, 
01280                        int NumImportIDs,
01281                        int * ImportLIDs, 
01282                        int LenImports,
01283                        char * Imports,
01284                        int & SizeOfPacket, 
01285                        Epetra_Distributor & Distor,
01286                        Epetra_CombineMode CombineMode,
01287                        const Epetra_OffsetIndex * Indexor);
01288 
01290   int SortEntries();
01291 
01293   bool Sorted() const {return(Graph_->Sorted());};
01294 
01296   int MergeRedundantEntries();
01297 
01299   bool NoRedundancies() const {return(Graph_->NoRedundancies());};
01300 
01301   bool StaticGraph() const {return(StaticGraph_);};
01302 
01303   int GeneratePointObjects() const;
01304   int BlockMap2PointMap(const Epetra_BlockMap & BlockMap, Epetra_Map * & PointMap) const;
01305   int UpdateOperatorXY(const Epetra_MultiVector& X, const Epetra_MultiVector& Y) const;
01306 
01307   Epetra_CrsGraph * Graph_;
01308   bool Allocated_;
01309   bool StaticGraph_;
01310   bool UseTranspose_;
01311   bool constructedWithFilledGraph_;
01312   bool matrixFillCompleteCalled_;
01313   bool StorageOptimized_;
01314 
01315   int NumMyBlockRows_;
01316 
01317   Epetra_DataAccess CV_;
01318 
01319 
01320   int * NumBlockEntriesPerRow_;
01321   int * NumAllocatedBlockEntriesPerRow_;
01322   int ** Indices_;
01323   int * ElementSizeList_;
01324   int * FirstPointInElementList_;
01325 
01326   Epetra_SerialDenseMatrix ***Entries_;
01327 
01328   double *All_Values_Orig_;
01329   double *All_Values_;
01330 
01331   mutable double NormInf_;
01332   mutable double NormOne_;
01333   mutable double NormFrob_;
01334 
01335   mutable Epetra_MultiVector * ImportVector_;
01336   mutable Epetra_MultiVector * ExportVector_;
01337 
01338   // State variables needed for constructing matrix entry-by-entry
01339   mutable int *TempRowDims_;
01340   mutable Epetra_SerialDenseMatrix **TempEntries_;
01341   mutable int LenTemps_;
01342   mutable int CurBlockRow_;
01343   mutable int CurNumBlockEntries_;
01344   mutable int * CurBlockIndices_;
01345   mutable int CurEntry_;
01346   mutable bool CurIndicesAreLocal_;
01347   mutable Epetra_CombineMode CurSubmitMode_;
01348   
01349   // State variables needed for extracting entries
01350   mutable int CurExtractBlockRow_;
01351   mutable int CurExtractEntry_; 
01352   mutable int CurExtractNumBlockEntries_;
01353   mutable bool CurExtractIndicesAreLocal_;
01354   mutable bool CurExtractView_;
01355   mutable int CurRowDim_;
01356 
01357   // State variable for extracting block diagonal entries
01358   mutable int CurBlockDiag_;
01359 
01360   // Maps and importer that support the Epetra_RowMatrix interface
01361   mutable Epetra_Map * RowMatrixRowMap_;
01362   mutable Epetra_Map * RowMatrixColMap_;
01363   mutable Epetra_Import * RowMatrixImporter_;
01364 
01365   // Maps that support the Epetra_Operator interface
01366   mutable Epetra_Map * OperatorDomainMap_;
01367   mutable Epetra_Map * OperatorRangeMap_;
01368   mutable Epetra_MultiVector * OperatorX_;
01369   mutable Epetra_MultiVector * OperatorY_;
01370 
01371   // bool to indicate if above four point maps and importer have already been created
01372   mutable bool HavePointObjects_;
01373 
01374   bool squareFillCompleteCalled_;
01375 };
01376 
01377 #endif /* EPETRA_VBRMATRIX_H */

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