EpetraExt_HypreIJMatrix.h

Go to the documentation of this file.

//@HEADER
// ***********************************************************************
//
//     EpetraExt: Epetra Extended - Linear Algebra Services Package
//                 Copyright (2011) Sandia Corporation
//
// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
// the U.S. Government retains certain rights in this software.
//
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are
// met:
//
// 1. Redistributions of source code must retain the above copyright
// notice, this list of conditions and the following disclaimer.
//
// 2. Redistributions in binary form must reproduce the above copyright
// notice, this list of conditions and the following disclaimer in the
// documentation and/or other materials provided with the distribution.
//
// 3. Neither the name of the Corporation nor the names of the
// contributors may be used to endorse or promote products derived from
// this software without specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
//
// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
//
// ***********************************************************************
//@HEADER

#ifndef EPETRAEXT_HYPREIJMATRIX_H_
#define EPETRAEXT_HYPREIJMATRIX_H_

// Trilinos source files
#include "Epetra_Object.h"
#include "Epetra_CompObject.h"
#include "Epetra_BasicRowMatrix.h"
#include "Epetra_Map.h"
#include "Epetra_Import.h"
#include "Epetra_MpiComm.h"

//Hypre source files
#include "krylov.h"
#include "HYPRE_parcsr_ls.h"
#include "_hypre_parcsr_mv.h"
#include "HYPRE_parcsr_mv.h"
#include "HYPRE_IJ_mv.h"
#include "_hypre_IJ_mv.h"
#include "HYPRE.h"

class Epetra_Vector;
class Epetra_MultiVector;
class Epetra_Import;


#ifndef HYPRE_ENUMS
#define HYPRE_ENUMS

  enum Hypre_Solver{ 
    BoomerAMG, 
    ParaSails, 
    Euclid,    
    AMS,       
    Hybrid,    
    PCG,       
    GMRES,     
    FlexGMRES, 
    LGMRES,    
    BiCGSTAB   
  };

  enum Hypre_Chooser{
    Solver,        
    Preconditioner 
  };
#endif //HYPRE_ENUMS

class EpetraExt_HypreIJMatrix: public Epetra_BasicRowMatrix  {
      
public:



  EpetraExt_HypreIJMatrix(HYPRE_IJMatrix matrix);

  virtual ~EpetraExt_HypreIJMatrix();
  



  int ExtractMyRowCopy(int MyRow, int Length, int & NumEntries, double *Values, int * Indices) const;
    

  int ExtractMyEntryView(int CurEntry, double *&Value, int &RowIndex, int &ColIndex);
    

  int ExtractMyEntryView(int CurEntry, const double *&Value, int &RowIndex, int &ColIndex) const;




  int SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, int), int parameter);


  int SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, double), double parameter);


  int SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, double, int), double parameter1, int parameter2);


  int SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, int, int), int parameter1, int parameter2);


  int SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, double*), double* parameter);


  int SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, int*), int* parameter);


  int SetParameter(Hypre_Chooser chooser, Hypre_Solver Solver, bool transpose=false);


  int SetParameter(bool UsePreconditioner);


  int SetParameter(Hypre_Chooser chooser);



  int Multiply(bool TransA, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;


  int Solve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;


  int LeftScale(const Epetra_Vector& X);



  int RightScale(const Epetra_Vector& X);




  int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const {
    return(EpetraExt_HypreIJMatrix::Multiply(EpetraExt_HypreIJMatrix::UseTranspose(), X, Y));};


  int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const {
    return(EpetraExt_HypreIJMatrix::Solve(EpetraExt_HypreIJMatrix::UpperTriangular(), EpetraExt_HypreIJMatrix::UseTranspose(), false, X, Y));};

  virtual bool UseTranspose() const {return(false);}




  int NumMyRowEntries(int MyRow, int & NumEntries) const;

  HYPRE_IJMatrix& GetMatrix(){ return Matrix_;}

protected:

  int InitializeDefaults();
   
  // These methods are needed only because the create methods in Hypre sometimes take an MPI_Comm but not always. 
  // They simply call the create solver in the correct way.
  int Hypre_BoomerAMGCreate(MPI_Comm comm, HYPRE_Solver *solver)
    { return HYPRE_BoomerAMGCreate(solver);}

  int Hypre_ParaSailsCreate(MPI_Comm comm, HYPRE_Solver *solver)
    { return HYPRE_ParaSailsCreate(comm, solver);}

  int Hypre_EuclidCreate(MPI_Comm comm, HYPRE_Solver *solver)
    { return HYPRE_EuclidCreate(comm, solver);}

  int Hypre_AMSCreate(MPI_Comm comm, HYPRE_Solver *solver)
    { return HYPRE_AMSCreate(solver);}

  int Hypre_ParCSRHybridCreate(MPI_Comm comm, HYPRE_Solver *solver)
    { return HYPRE_ParCSRHybridCreate(solver);}

  int Hypre_ParCSRPCGCreate(MPI_Comm comm, HYPRE_Solver *solver)
    { return HYPRE_ParCSRPCGCreate(comm, solver);}

  int Hypre_ParCSRGMRESCreate(MPI_Comm comm, HYPRE_Solver *solver)
    { return HYPRE_ParCSRGMRESCreate(comm, solver);}

  int Hypre_ParCSRFlexGMRESCreate(MPI_Comm comm, HYPRE_Solver *solver)
    { return HYPRE_ParCSRFlexGMRESCreate(comm, solver);}

  int Hypre_ParCSRLGMRESCreate(MPI_Comm comm, HYPRE_Solver *solver)
    { return HYPRE_ParCSRLGMRESCreate(comm, solver);}

  int Hypre_ParCSRBiCGSTABCreate(MPI_Comm comm, HYPRE_Solver *solver)
    { return HYPRE_ParCSRBiCGSTABCreate(comm, solver);}

  int CreateSolver();
  int CreatePrecond();
  // These two methods setup the solver or preconditioner by calling the pointer. 
  // They are const because they are called from the Solve() routine.
  // They really aren't const because they change the value of IsSolverSetup_. This is because it should only be called if it isn't setup.
  int SetupSolver() const;
  int SetupPrecond() const;

  // Hypre variables
  mutable HYPRE_IJMatrix Matrix_;
  mutable HYPRE_ParCSRMatrix ParMatrix_;
  mutable HYPRE_IJVector X_hypre;
  mutable HYPRE_IJVector Y_hypre;
  mutable HYPRE_ParVector par_x;
  mutable HYPRE_ParVector par_y;
  mutable hypre_ParVector *x_vec;
  mutable hypre_ParVector *y_vec;
  mutable hypre_Vector *x_local;
  mutable hypre_Vector *y_local;
  mutable HYPRE_Solver Solver_;
  mutable HYPRE_Solver Preconditioner_;
  // The following are pointers to functions to use the solver and preconditioner.
  int (EpetraExt_HypreIJMatrix::*SolverCreatePtr_)(MPI_Comm, HYPRE_Solver*);
  int (*SolverDestroyPtr_)(HYPRE_Solver);
  int (*SolverSetupPtr_)(HYPRE_Solver, HYPRE_ParCSRMatrix, HYPRE_ParVector, HYPRE_ParVector);
  int (*SolverSolvePtr_)(HYPRE_Solver, HYPRE_ParCSRMatrix, HYPRE_ParVector, HYPRE_ParVector);
  int (*SolverPrecondPtr_)(HYPRE_Solver, HYPRE_PtrToParSolverFcn, HYPRE_PtrToParSolverFcn, HYPRE_Solver);
  int (EpetraExt_HypreIJMatrix::*PrecondCreatePtr_)(MPI_Comm, HYPRE_Solver*);
  int (*PrecondDestroyPtr_)(HYPRE_Solver);
  int (*PrecondSetupPtr_)(HYPRE_Solver, HYPRE_ParCSRMatrix, HYPRE_ParVector, HYPRE_ParVector);
  int (*PrecondSolvePtr_)(HYPRE_Solver, HYPRE_ParCSRMatrix, HYPRE_ParVector, HYPRE_ParVector);
    
  int NumMyRows_;
  int NumGlobalRows_;
  int NumGlobalCols_;
  int MyRowStart_;
  int MyRowEnd_;
  mutable int MatType_;
  bool TransposeSolve_;
  Hypre_Chooser SolveOrPrec_; 
  bool *IsSolverSetup_; 
  bool *IsPrecondSetup_;
};
#endif /* EPETRAEXT_HYPREIJMATRIX_H_ */