IFPACK  Development
Ifpack_Hypre.cpp
1 /*@HEADER
2 // ***********************************************************************
3 //
4 // Ifpack: Object-Oriented Algebraic Preconditioner Package
5 // Copyright (2002) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 //@HEADER
41 */
42 #include "Ifpack_Hypre.h"
43 #if defined(HAVE_HYPRE) && defined(HAVE_MPI)
44 
45 #include "Ifpack_Utils.h"
46 #include "Epetra_MpiComm.h"
47 #include "Epetra_IntVector.h"
48 #include "Epetra_Import.h"
49 #include "Teuchos_ParameterList.hpp"
50 #include "Teuchos_RCP.hpp"
51 
52 using Teuchos::RCP;
53 using Teuchos::rcp;
54 
55 Ifpack_Hypre::Ifpack_Hypre(Epetra_RowMatrix* A):
56  A_(rcp(A,false)),
57  UseTranspose_(false),
58  IsInitialized_(false),
59  IsComputed_(false),
60  Label_(),
61  NumInitialize_(0),
62  NumCompute_(0),
63  NumApplyInverse_(0),
64  InitializeTime_(0.0),
65  ComputeTime_(0.0),
66  ApplyInverseTime_(0.0),
67  ComputeFlops_(0.0),
68  ApplyInverseFlops_(0.0),
69  Time_(A_->Comm()),
70  SolveOrPrec_(Solver),
71  NumFunsToCall_(0),
72  SolverType_(PCG),
73  PrecondType_(Euclid),
74  UsePreconditioner_(false),
75  NiceRowMap_(true)
76 {
77  IsSolverSetup_ = new bool[1];
78  IsPrecondSetup_ = new bool[1];
79  IsSolverSetup_[0] = false;
80  IsPrecondSetup_[0] = false;
81  MPI_Comm comm = GetMpiComm();
82  int ilower = A_->RowMatrixRowMap().MinMyGID();
83  int iupper = A_->RowMatrixRowMap().MaxMyGID();
84  // Need to check if the RowMap is the way Hypre expects (if not more difficult)
85  std::vector<int> ilowers; ilowers.resize(Comm().NumProc());
86  std::vector<int> iuppers; iuppers.resize(Comm().NumProc());
87  int myLower[1]; myLower[0] = ilower;
88  int myUpper[1]; myUpper[0] = iupper;
89  Comm().GatherAll(myLower, &ilowers[0], 1);
90  Comm().GatherAll(myUpper, &iuppers[0], 1);
91  for(int i = 0; i < Comm().NumProc()-1; i++){
92  NiceRowMap_ = (NiceRowMap_ && iuppers[i]+1 == ilowers[i+1]);
93  }
94  if(!NiceRowMap_){
95  ilower = (A_->NumGlobalRows() / Comm().NumProc())*Comm().MyPID();
96  iupper = (A_->NumGlobalRows() / Comm().NumProc())*(Comm().MyPID()+1)-1;
97  if(Comm().MyPID() == Comm().NumProc()-1){
98  iupper = A_-> NumGlobalRows()-1;
99  }
100  }
101 
102  // Next create vectors that will be used when ApplyInverse() is called
103  IFPACK_CHK_ERRV(HYPRE_IJVectorCreate(comm, ilower, iupper, &XHypre_));
104  IFPACK_CHK_ERRV(HYPRE_IJVectorSetObjectType(XHypre_, HYPRE_PARCSR));
105  IFPACK_CHK_ERRV(HYPRE_IJVectorInitialize(XHypre_));
106  IFPACK_CHK_ERRV(HYPRE_IJVectorAssemble(XHypre_));
107  IFPACK_CHK_ERRV(HYPRE_IJVectorGetObject(XHypre_, (void**) &ParX_));
108 
109  IFPACK_CHK_ERRV(HYPRE_IJVectorCreate(comm, ilower, iupper, &YHypre_));
110  IFPACK_CHK_ERRV(HYPRE_IJVectorSetObjectType(YHypre_, HYPRE_PARCSR));
111  IFPACK_CHK_ERRV(HYPRE_IJVectorInitialize(YHypre_));
112  IFPACK_CHK_ERRV(HYPRE_IJVectorAssemble(YHypre_));
113  IFPACK_CHK_ERRV(HYPRE_IJVectorGetObject(YHypre_, (void**) &ParY_));
114 
115  XVec_ = (hypre_ParVector *) hypre_IJVectorObject(((hypre_IJVector *) XHypre_));
116  XLocal_ = hypre_ParVectorLocalVector(XVec_);
117 
118  YVec_ = (hypre_ParVector *) hypre_IJVectorObject(((hypre_IJVector *) YHypre_));
119  YLocal_ = hypre_ParVectorLocalVector(YVec_);
120  std::vector<int> rows; rows.resize(iupper - ilower +1);
121  for(int i = ilower; i <= iupper; i++){
122  rows[i-ilower] = i;
123  }
124  MySimpleMap_ = rcp(new Epetra_Map(-1, iupper-ilower+1, &rows[0], 0, Comm()));
125 } //Constructor
126 
127 //==============================================================================
128 void Ifpack_Hypre::Destroy(){
129  if(IsInitialized()){
130  IFPACK_CHK_ERRV(HYPRE_IJMatrixDestroy(HypreA_));
131  }
132  IFPACK_CHK_ERRV(HYPRE_IJVectorDestroy(XHypre_));
133  IFPACK_CHK_ERRV(HYPRE_IJVectorDestroy(YHypre_));
134  if(IsSolverSetup_[0]){
135  IFPACK_CHK_ERRV(SolverDestroyPtr_(Solver_));
136  }
137  if(IsPrecondSetup_[0]){
138  IFPACK_CHK_ERRV(PrecondDestroyPtr_(Preconditioner_));
139  }
140  delete[] IsSolverSetup_;
141  delete[] IsPrecondSetup_;
142 } //Destroy()
143 
144 //==============================================================================
145 int Ifpack_Hypre::Initialize(){
146  Time_.ResetStartTime();
147  MPI_Comm comm = GetMpiComm();
148  int ilower = MySimpleMap_->MinMyGID();
149  int iupper = MySimpleMap_->MaxMyGID();
150  IFPACK_CHK_ERR(HYPRE_IJMatrixCreate(comm, ilower, iupper, ilower, iupper, &HypreA_));
151  IFPACK_CHK_ERR(HYPRE_IJMatrixSetObjectType(HypreA_, HYPRE_PARCSR));
152  IFPACK_CHK_ERR(HYPRE_IJMatrixInitialize(HypreA_));
153  for(int i = 0; i < A_->NumMyRows(); i++){
154  int numElements;
155  IFPACK_CHK_ERR(A_->NumMyRowEntries(i,numElements));
156  std::vector<int> indices; indices.resize(numElements);
157  std::vector<double> values; values.resize(numElements);
158  int numEntries;
159  IFPACK_CHK_ERR(A_->ExtractMyRowCopy(i, numElements, numEntries, &values[0], &indices[0]));
160  for(int j = 0; j < numEntries; j++){
161  indices[j] = A_->RowMatrixColMap().GID(indices[j]);
162  }
163  int GlobalRow[1];
164  GlobalRow[0] = A_->RowMatrixRowMap().GID(i);
165  IFPACK_CHK_ERR(HYPRE_IJMatrixSetValues(HypreA_, 1, &numEntries, GlobalRow, &indices[0], &values[0]));
166  }
167  IFPACK_CHK_ERR(HYPRE_IJMatrixAssemble(HypreA_));
168  IFPACK_CHK_ERR(HYPRE_IJMatrixGetObject(HypreA_, (void**)&ParMatrix_));
169  IsInitialized_=true;
170  NumInitialize_ = NumInitialize_ + 1;
171  InitializeTime_ = InitializeTime_ + Time_.ElapsedTime();
172  return 0;
173 } //Initialize()
174 
175 //==============================================================================
176 int Ifpack_Hypre::SetParameters(Teuchos::ParameterList& list){
177  List_ = list;
178  Hypre_Solver solType = list.get("Solver", PCG);
179  SolverType_ = solType;
180  Hypre_Solver precType = list.get("Preconditioner", Euclid);
181  PrecondType_ = precType;
182  Hypre_Chooser chooser = list.get("SolveOrPrecondition", Solver);
183  SolveOrPrec_ = chooser;
184  bool SetPrecond = list.get("SetPreconditioner", false);
185  IFPACK_CHK_ERR(SetParameter(SetPrecond));
186  int NumFunctions = list.get("NumFunctions", 0);
187  FunsToCall_.clear();
188  NumFunsToCall_ = 0;
189  if(NumFunctions > 0){
190  RCP<FunctionParameter>* params = list.get<RCP<FunctionParameter>*>("Functions");
191  for(int i = 0; i < NumFunctions; i++){
192  IFPACK_CHK_ERR(AddFunToList(params[i]));
193  }
194  }
195  return 0;
196 } //SetParameters()
197 
198 //==============================================================================
199 int Ifpack_Hypre::AddFunToList(RCP<FunctionParameter> NewFun){
200  NumFunsToCall_ = NumFunsToCall_+1;
201  FunsToCall_.resize(NumFunsToCall_);
202  FunsToCall_[NumFunsToCall_-1] = NewFun;
203  return 0;
204 } //AddFunToList()
205 
206 //==============================================================================
207 int Ifpack_Hypre::SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, int), int parameter){
208  RCP<FunctionParameter> temp = rcp(new FunctionParameter(chooser, pt2Func, parameter));
209  IFPACK_CHK_ERR(AddFunToList(temp));
210  return 0;
211 } //SetParameter() - int function pointer
212 
213 //==============================================================================
214 int Ifpack_Hypre::SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, double), double parameter){
215  RCP<FunctionParameter> temp = rcp(new FunctionParameter(chooser, pt2Func, parameter));
216  IFPACK_CHK_ERR(AddFunToList(temp));
217  return 0;
218 } //SetParameter() - double function pointer
219 
220 //==============================================================================
221 int Ifpack_Hypre::SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, double, int), double parameter1, int parameter2){
222  RCP<FunctionParameter> temp = rcp(new FunctionParameter(chooser, pt2Func, parameter1, parameter2));
223  IFPACK_CHK_ERR(AddFunToList(temp));
224  return 0;
225 } //SetParameter() - double,int function pointer
226 
227 //==============================================================================
228 int Ifpack_Hypre::SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, int, int), int parameter1, int parameter2){
229  RCP<FunctionParameter> temp = rcp(new FunctionParameter(chooser, pt2Func, parameter1, parameter2));
230  IFPACK_CHK_ERR(AddFunToList(temp));
231  return 0;
232 } //SetParameter() int,int function pointer
233 
234 //==============================================================================
235 int Ifpack_Hypre::SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, double*), double* parameter){
236  RCP<FunctionParameter> temp = rcp(new FunctionParameter(chooser, pt2Func, parameter));
237  IFPACK_CHK_ERR(AddFunToList(temp));
238  return 0;
239 } //SetParameter() - double* function pointer
240 
241 //==============================================================================
242 int Ifpack_Hypre::SetParameter(Hypre_Chooser chooser, int (*pt2Func)(HYPRE_Solver, int*), int* parameter){
243  RCP<FunctionParameter> temp = rcp(new FunctionParameter(chooser, pt2Func, parameter));
244  IFPACK_CHK_ERR(AddFunToList(temp));
245  return 0;
246 } //SetParameter() - int* function pointer
247 
248 //==============================================================================
249 int Ifpack_Hypre::SetParameter(Hypre_Chooser chooser, Hypre_Solver solver){
250  if(chooser == Solver){
251  SolverType_ = solver;
252  } else {
253  PrecondType_ = solver;
254  }
255  return 0;
256 } //SetParameter() - set type of solver
257 
258 //==============================================================================
259 int Ifpack_Hypre::Compute(){
260  if(IsInitialized() == false){
261  IFPACK_CHK_ERR(Initialize());
262  }
263  Time_.ResetStartTime();
264  IFPACK_CHK_ERR(SetSolverType(SolverType_));
265  IFPACK_CHK_ERR(SetPrecondType(PrecondType_));
266  CallFunctions();
267  if(UsePreconditioner_){
268  if(SolverPrecondPtr_ != NULL){
269  IFPACK_CHK_ERR(SolverPrecondPtr_(Solver_, PrecondSolvePtr_, PrecondSetupPtr_, Preconditioner_));
270  }
271  }
272  if(SolveOrPrec_ == Solver){
273  IFPACK_CHK_ERR(SolverSetupPtr_(Solver_, ParMatrix_, ParX_, ParY_));
274  IsSolverSetup_[0] = true;
275  } else {
276  IFPACK_CHK_ERR(PrecondSetupPtr_(Preconditioner_, ParMatrix_, ParX_, ParY_));
277  IsPrecondSetup_[0] = true;
278  }
279  IsComputed_ = true;
280  NumCompute_ = NumCompute_ + 1;
281  ComputeTime_ = ComputeTime_ + Time_.ElapsedTime();
282  return 0;
283 } //Compute()
284 
285 //==============================================================================
286 int Ifpack_Hypre::CallFunctions() const{
287  for(int i = 0; i < NumFunsToCall_; i++){
288  IFPACK_CHK_ERR(FunsToCall_[i]->CallFunction(Solver_, Preconditioner_));
289  }
290  return 0;
291 } //CallFunctions()
292 
293 //==============================================================================
294 int Ifpack_Hypre::ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const{
295  if(IsComputed() == false){
296  IFPACK_CHK_ERR(-1);
297  }
298  Time_.ResetStartTime();
299  bool SameVectors = false;
300  int NumVectors = X.NumVectors();
301  if (NumVectors != Y.NumVectors()) IFPACK_CHK_ERR(-1); // X and Y must have same number of vectors
302  if(X.Pointers() == Y.Pointers()){
303  SameVectors = true;
304  }
305  for(int VecNum = 0; VecNum < NumVectors; VecNum++) {
306  //Get values for current vector in multivector.
307  double * XValues;
308  IFPACK_CHK_ERR((*X(VecNum)).ExtractView(&XValues));
309  double * YValues;
310  if(!SameVectors){
311  IFPACK_CHK_ERR((*Y(VecNum)).ExtractView(&YValues));
312  } else {
313  YValues = new double[X.MyLength()];
314  }
315  // Temporarily make a pointer to data in Hypre for end
316  double *XTemp = XLocal_->data;
317  // Replace data in Hypre vectors with epetra values
318  XLocal_->data = XValues;
319  double *YTemp = YLocal_->data;
320  YLocal_->data = YValues;
321 
322  IFPACK_CHK_ERR(HYPRE_ParVectorSetConstantValues(ParY_, 0.0));
323  if(SolveOrPrec_ == Solver){
324  // Use the solver methods
325  IFPACK_CHK_ERR(SolverSolvePtr_(Solver_, ParMatrix_, ParX_, ParY_));
326  } else {
327  // Apply the preconditioner
328  IFPACK_CHK_ERR(PrecondSolvePtr_(Preconditioner_, ParMatrix_, ParX_, ParY_));
329  }
330  if(SameVectors){
331  int NumEntries = Y.MyLength();
332  std::vector<double> new_values; new_values.resize(NumEntries);
333  std::vector<int> new_indices; new_indices.resize(NumEntries);
334  for(int i = 0; i < NumEntries; i++){
335  new_values[i] = YValues[i];
336  new_indices[i] = i;
337  }
338  IFPACK_CHK_ERR((*Y(VecNum)).ReplaceMyValues(NumEntries, &new_values[0], &new_indices[0]));
339  delete[] YValues;
340  }
341  XLocal_->data = XTemp;
342  YLocal_->data = YTemp;
343  }
344  NumApplyInverse_ = NumApplyInverse_ + 1;
345  ApplyInverseTime_ = ApplyInverseTime_ + Time_.ElapsedTime();
346  return 0;
347 } //ApplyInverse()
348 
349 //==============================================================================
350 int Ifpack_Hypre::Multiply(bool TransA, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const{
351  if(IsInitialized() == false){
352  IFPACK_CHK_ERR(-1);
353  }
354  bool SameVectors = false;
355  int NumVectors = X.NumVectors();
356  if (NumVectors != Y.NumVectors()) IFPACK_CHK_ERR(-1); // X and Y must have same number of vectors
357  if(X.Pointers() == Y.Pointers()){
358  SameVectors = true;
359  }
360  for(int VecNum = 0; VecNum < NumVectors; VecNum++) {
361  //Get values for current vector in multivector.
362  double * XValues;
363  double * YValues;
364  IFPACK_CHK_ERR((*X(VecNum)).ExtractView(&XValues));
365  double *XTemp = XLocal_->data;
366  double *YTemp = YLocal_->data;
367  if(!SameVectors){
368  IFPACK_CHK_ERR((*Y(VecNum)).ExtractView(&YValues));
369  } else {
370  YValues = new double[X.MyLength()];
371  }
372  YLocal_->data = YValues;
373  IFPACK_CHK_ERR(HYPRE_ParVectorSetConstantValues(ParY_,0.0));
374  // Temporarily make a pointer to data in Hypre for end
375  // Replace data in Hypre vectors with epetra values
376  XLocal_->data = XValues;
377  // Do actual computation.
378  if(TransA) {
379  // Use transpose of A in multiply
380  IFPACK_CHK_ERR(HYPRE_ParCSRMatrixMatvecT(1.0, ParMatrix_, ParX_, 1.0, ParY_));
381  } else {
382  IFPACK_CHK_ERR(HYPRE_ParCSRMatrixMatvec(1.0, ParMatrix_, ParX_, 1.0, ParY_));
383  }
384  if(SameVectors){
385  int NumEntries = Y.MyLength();
386  std::vector<double> new_values; new_values.resize(NumEntries);
387  std::vector<int> new_indices; new_indices.resize(NumEntries);
388  for(int i = 0; i < NumEntries; i++){
389  new_values[i] = YValues[i];
390  new_indices[i] = i;
391  }
392  IFPACK_CHK_ERR((*Y(VecNum)).ReplaceMyValues(NumEntries, &new_values[0], &new_indices[0]));
393  delete[] YValues;
394  }
395  XLocal_->data = XTemp;
396  YLocal_->data = YTemp;
397  }
398  return 0;
399 } //Multiply()
400 
401 //==============================================================================
402 std::ostream& Ifpack_Hypre::Print(std::ostream& os) const{
403  using std::endl;
404 
405  if (!Comm().MyPID()) {
406  os << endl;
407  os << "================================================================================" << endl;
408  os << "Ifpack_Hypre: " << Label () << endl << endl;
409  os << "Using " << Comm().NumProc() << " processors." << endl;
410  os << "Global number of rows = " << A_->NumGlobalRows() << endl;
411  os << "Global number of nonzeros = " << A_->NumGlobalNonzeros() << endl;
412  os << "Condition number estimate = " << Condest() << endl;
413  os << endl;
414  os << "Phase # calls Total Time (s) Total MFlops MFlops/s" << endl;
415  os << "----- ------- -------------- ------------ --------" << endl;
416  os << "Initialize() " << std::setw(5) << NumInitialize_
417  << " " << std::setw(15) << InitializeTime_
418  << " 0.0 0.0" << endl;
419  os << "Compute() " << std::setw(5) << NumCompute_
420  << " " << std::setw(15) << ComputeTime_
421  << " " << std::setw(15) << 1.0e-6 * ComputeFlops_;
422  if (ComputeTime_ != 0.0)
423  os << " " << std::setw(15) << 1.0e-6 * ComputeFlops_ / ComputeTime_ << endl;
424  else
425  os << " " << std::setw(15) << 0.0 << endl;
426  os << "ApplyInverse() " << std::setw(5) << NumApplyInverse_
427  << " " << std::setw(15) << ApplyInverseTime_
428  << " " << std::setw(15) << 1.0e-6 * ApplyInverseFlops_;
429  if (ApplyInverseTime_ != 0.0)
430  os << " " << std::setw(15) << 1.0e-6 * ApplyInverseFlops_ / ApplyInverseTime_ << endl;
431  else
432  os << " " << std::setw(15) << 0.0 << endl;
433  os << "================================================================================" << endl;
434  os << endl;
435  }
436  return os;
437 } //Print()
438 
439 //==============================================================================
440 double Ifpack_Hypre::Condest(const Ifpack_CondestType CT,
441  const int MaxIters,
442  const double Tol,
443  Epetra_RowMatrix* Matrix_in){
444  if (!IsComputed()) // cannot compute right now
445  return(-1.0);
446  Condest_ = Ifpack_Condest(*this, CT, MaxIters, Tol, Matrix_in);
447  return(Condest_);
448 } //Condest()
449 
450 //==============================================================================
451 int Ifpack_Hypre::SetSolverType(Hypre_Solver Solver){
452  switch(Solver) {
453  case BoomerAMG:
454  if(IsSolverSetup_[0]){
455  SolverDestroyPtr_(Solver_);
456  IsSolverSetup_[0] = false;
457  }
458  SolverCreatePtr_ = &Ifpack_Hypre::Hypre_BoomerAMGCreate;
459  SolverDestroyPtr_ = &HYPRE_BoomerAMGDestroy;
460  SolverSetupPtr_ = &HYPRE_BoomerAMGSetup;
461  SolverPrecondPtr_ = NULL;
462  SolverSolvePtr_ = &HYPRE_BoomerAMGSolve;
463  break;
464  case AMS:
465  if(IsSolverSetup_[0]){
466  SolverDestroyPtr_(Solver_);
467  IsSolverSetup_[0] = false;
468  }
469  SolverCreatePtr_ = &Ifpack_Hypre::Hypre_AMSCreate;
470  SolverDestroyPtr_ = &HYPRE_AMSDestroy;
471  SolverSetupPtr_ = &HYPRE_AMSSetup;
472  SolverSolvePtr_ = &HYPRE_AMSSolve;
473  SolverPrecondPtr_ = NULL;
474  break;
475  case Hybrid:
476  if(IsSolverSetup_[0]){
477  SolverDestroyPtr_(Solver_);
478  IsSolverSetup_[0] = false;
479  }
480  SolverCreatePtr_ = &Ifpack_Hypre::Hypre_ParCSRHybridCreate;
481  SolverDestroyPtr_ = &HYPRE_ParCSRHybridDestroy;
482  SolverSetupPtr_ = &HYPRE_ParCSRHybridSetup;
483  SolverSolvePtr_ = &HYPRE_ParCSRHybridSolve;
484  SolverPrecondPtr_ = &HYPRE_ParCSRHybridSetPrecond;
485  break;
486  case PCG:
487  if(IsSolverSetup_[0]){
488  SolverDestroyPtr_(Solver_);
489  IsSolverSetup_[0] = false;
490  }
491  SolverCreatePtr_ = &Ifpack_Hypre::Hypre_ParCSRPCGCreate;
492  SolverDestroyPtr_ = &HYPRE_ParCSRPCGDestroy;
493  SolverSetupPtr_ = &HYPRE_ParCSRPCGSetup;
494  SolverSolvePtr_ = &HYPRE_ParCSRPCGSolve;
495  SolverPrecondPtr_ = &HYPRE_ParCSRPCGSetPrecond;
496  break;
497  case GMRES:
498  if(IsSolverSetup_[0]){
499  SolverDestroyPtr_(Solver_);
500  IsSolverSetup_[0] = false;
501  }
502  SolverCreatePtr_ = &Ifpack_Hypre::Hypre_ParCSRGMRESCreate;
503  SolverDestroyPtr_ = &HYPRE_ParCSRGMRESDestroy;
504  SolverSetupPtr_ = &HYPRE_ParCSRGMRESSetup;
505  SolverPrecondPtr_ = &HYPRE_ParCSRGMRESSetPrecond;
506  break;
507  case FlexGMRES:
508  if(IsSolverSetup_[0]){
509  SolverDestroyPtr_(Solver_);
510  IsSolverSetup_[0] = false;
511  }
512  SolverCreatePtr_ = &Ifpack_Hypre::Hypre_ParCSRFlexGMRESCreate;
513  SolverDestroyPtr_ = &HYPRE_ParCSRFlexGMRESDestroy;
514  SolverSetupPtr_ = &HYPRE_ParCSRFlexGMRESSetup;
515  SolverSolvePtr_ = &HYPRE_ParCSRFlexGMRESSolve;
516  SolverPrecondPtr_ = &HYPRE_ParCSRFlexGMRESSetPrecond;
517  break;
518  case LGMRES:
519  if(IsSolverSetup_[0]){
520  SolverDestroyPtr_(Solver_);
521  IsSolverSetup_[0] = false;
522  }
523  SolverCreatePtr_ = &Ifpack_Hypre::Hypre_ParCSRLGMRESCreate;
524  SolverDestroyPtr_ = &HYPRE_ParCSRLGMRESDestroy;
525  SolverSetupPtr_ = &HYPRE_ParCSRLGMRESSetup;
526  SolverSolvePtr_ = &HYPRE_ParCSRLGMRESSolve;
527  SolverPrecondPtr_ = &HYPRE_ParCSRLGMRESSetPrecond;
528  break;
529  case BiCGSTAB:
530  if(IsSolverSetup_[0]){
531  SolverDestroyPtr_(Solver_);
532  IsSolverSetup_[0] = false;
533  }
534  SolverCreatePtr_ = &Ifpack_Hypre::Hypre_ParCSRBiCGSTABCreate;
535  SolverDestroyPtr_ = &HYPRE_ParCSRBiCGSTABDestroy;
536  SolverSetupPtr_ = &HYPRE_ParCSRBiCGSTABSetup;
537  SolverSolvePtr_ = &HYPRE_ParCSRBiCGSTABSolve;
538  SolverPrecondPtr_ = &HYPRE_ParCSRBiCGSTABSetPrecond;
539  break;
540  default:
541  return -1;
542  }
543  CreateSolver();
544  return 0;
545 } //SetSolverType()
546 
547 //==============================================================================
548 int Ifpack_Hypre::SetPrecondType(Hypre_Solver Precond){
549  switch(Precond) {
550  case BoomerAMG:
551  if(IsPrecondSetup_[0]){
552  PrecondDestroyPtr_(Preconditioner_);
553  IsPrecondSetup_[0] = false;
554  }
555  PrecondCreatePtr_ = &Ifpack_Hypre::Hypre_BoomerAMGCreate;
556  PrecondDestroyPtr_ = &HYPRE_BoomerAMGDestroy;
557  PrecondSetupPtr_ = &HYPRE_BoomerAMGSetup;
558  PrecondSolvePtr_ = &HYPRE_BoomerAMGSolve;
559  break;
560  case ParaSails:
561  if(IsPrecondSetup_[0]){
562  PrecondDestroyPtr_(Preconditioner_);
563  IsPrecondSetup_[0] = false;
564  }
565  PrecondCreatePtr_ = &Ifpack_Hypre::Hypre_ParaSailsCreate;
566  PrecondDestroyPtr_ = &HYPRE_ParaSailsDestroy;
567  PrecondSetupPtr_ = &HYPRE_ParaSailsSetup;
568  PrecondSolvePtr_ = &HYPRE_ParaSailsSolve;
569  break;
570  case Euclid:
571  if(IsPrecondSetup_[0]){
572  PrecondDestroyPtr_(Preconditioner_);
573  IsPrecondSetup_[0] = false;
574  }
575  PrecondCreatePtr_ = &Ifpack_Hypre::Hypre_EuclidCreate;
576  PrecondDestroyPtr_ = &HYPRE_EuclidDestroy;
577  PrecondSetupPtr_ = &HYPRE_EuclidSetup;
578  PrecondSolvePtr_ = &HYPRE_EuclidSolve;
579  break;
580  case AMS:
581  if(IsPrecondSetup_[0]){
582  PrecondDestroyPtr_(Preconditioner_);
583  IsPrecondSetup_[0] = false;
584  }
585  PrecondCreatePtr_ = &Ifpack_Hypre::Hypre_AMSCreate;
586  PrecondDestroyPtr_ = &HYPRE_AMSDestroy;
587  PrecondSetupPtr_ = &HYPRE_AMSSetup;
588  PrecondSolvePtr_ = &HYPRE_AMSSolve;
589  break;
590  default:
591  return -1;
592  }
593  CreatePrecond();
594  return 0;
595 
596 } //SetPrecondType()
597 
598 //==============================================================================
599 int Ifpack_Hypre::CreateSolver(){
600  MPI_Comm comm;
601  HYPRE_ParCSRMatrixGetComm(ParMatrix_, &comm);
602  return (this->*SolverCreatePtr_)(comm, &Solver_);
603 } //CreateSolver()
604 
605 //==============================================================================
606 int Ifpack_Hypre::CreatePrecond(){
607  MPI_Comm comm;
608  HYPRE_ParCSRMatrixGetComm(ParMatrix_, &comm);
609  return (this->*PrecondCreatePtr_)(comm, &Preconditioner_);
610 } //CreatePrecond()
611 
612 #endif // HAVE_HYPRE && HAVE_MPI