1 2 /* 3 Provides an interface to the LLNL package hypre 4 */ 5 6 /* Must use hypre 2.0.0 or more recent. */ 7 8 #include <petsc/private/pcimpl.h> /*I "petscpc.h" I*/ 9 #include <../src/dm/impls/da/hypre/mhyp.h> 10 #include <_hypre_parcsr_ls.h> 11 12 static PetscBool cite = PETSC_FALSE; 13 static const char hypreCitation[] = "@manual{hypre-web-page,\n title = {{\\sl hypre}: High Performance Preconditioners},\n organization = {Lawrence Livermore National Laboratory},\n note = {\\url{http://www.llnl.gov/CASC/hypre/}}\n}\n"; 14 15 /* 16 Private context (data structure) for the preconditioner. 17 */ 18 typedef struct { 19 HYPRE_Solver hsolver; 20 HYPRE_IJMatrix ij; 21 HYPRE_IJVector b,x; 22 23 HYPRE_Int (*destroy)(HYPRE_Solver); 24 HYPRE_Int (*solve)(HYPRE_Solver,HYPRE_ParCSRMatrix,HYPRE_ParVector,HYPRE_ParVector); 25 HYPRE_Int (*setup)(HYPRE_Solver,HYPRE_ParCSRMatrix,HYPRE_ParVector,HYPRE_ParVector); 26 HYPRE_Int (*setdgrad)(HYPRE_Solver,HYPRE_ParCSRMatrix); 27 HYPRE_Int (*setdcurl)(HYPRE_Solver,HYPRE_ParCSRMatrix); 28 HYPRE_Int (*setcoord)(HYPRE_Solver,HYPRE_ParVector,HYPRE_ParVector,HYPRE_ParVector); 29 HYPRE_Int (*setdim)(HYPRE_Solver,HYPRE_Int); 30 31 MPI_Comm comm_hypre; 32 char *hypre_type; 33 34 /* options for Pilut and BoomerAMG*/ 35 PetscInt maxiter; 36 double tol; 37 38 /* options for Pilut */ 39 PetscInt factorrowsize; 40 41 /* options for ParaSails */ 42 PetscInt nlevels; 43 double threshhold; 44 double filter; 45 PetscInt sym; 46 double loadbal; 47 PetscInt logging; 48 PetscInt ruse; 49 PetscInt symt; 50 51 /* options for BoomerAMG */ 52 PetscBool printstatistics; 53 54 /* options for BoomerAMG */ 55 PetscInt cycletype; 56 PetscInt maxlevels; 57 double strongthreshold; 58 double maxrowsum; 59 PetscInt gridsweeps[3]; 60 PetscInt coarsentype; 61 PetscInt measuretype; 62 PetscInt relaxtype[3]; 63 double relaxweight; 64 double outerrelaxweight; 65 PetscInt relaxorder; 66 double truncfactor; 67 PetscBool applyrichardson; 68 PetscInt pmax; 69 PetscInt interptype; 70 PetscInt agg_nl; 71 PetscInt agg_num_paths; 72 PetscInt nodal_coarsen; 73 PetscBool nodal_relax; 74 PetscInt nodal_relax_levels; 75 76 /* options for AS (Auxiliary Space preconditioners) */ 77 PetscInt as_print; 78 PetscInt as_max_iter; 79 PetscReal as_tol; 80 PetscInt as_relax_type; 81 PetscInt as_relax_times; 82 PetscReal as_relax_weight; 83 PetscReal as_omega; 84 PetscInt as_amg_alpha_opts[5]; /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for vector Poisson (AMS) or Curl problem (ADS) */ 85 PetscReal as_amg_alpha_theta; /* AMG strength for vector Poisson (AMS) or Curl problem (ADS) */ 86 PetscInt as_amg_beta_opts[5]; /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for scalar Poisson (AMS) or vector Poisson (ADS) */ 87 PetscReal as_amg_beta_theta; /* AMG strength for scalar Poisson (AMS) or vector Poisson (ADS) */ 88 PetscInt ams_cycle_type; 89 PetscInt ads_cycle_type; 90 91 /* additional data */ 92 HYPRE_IJVector coords[3]; 93 HYPRE_IJVector constants[3]; 94 HYPRE_IJMatrix G; 95 HYPRE_IJMatrix C; 96 HYPRE_IJMatrix alpha_Poisson; 97 HYPRE_IJMatrix beta_Poisson; 98 PetscBool ams_beta_is_zero; 99 } PC_HYPRE; 100 101 #undef __FUNCT__ 102 #define __FUNCT__ "PCHYPREGetSolver" 103 PetscErrorCode PCHYPREGetSolver(PC pc,HYPRE_Solver *hsolver) 104 { 105 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 106 107 PetscFunctionBegin; 108 *hsolver = jac->hsolver; 109 PetscFunctionReturn(0); 110 } 111 112 #undef __FUNCT__ 113 #define __FUNCT__ "PCSetUp_HYPRE" 114 static PetscErrorCode PCSetUp_HYPRE(PC pc) 115 { 116 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 117 PetscErrorCode ierr; 118 HYPRE_ParCSRMatrix hmat; 119 HYPRE_ParVector bv,xv; 120 PetscInt bs; 121 122 PetscFunctionBegin; 123 if (!jac->hypre_type) { 124 ierr = PCHYPRESetType(pc,"boomeramg");CHKERRQ(ierr); 125 } 126 127 if (pc->setupcalled) { 128 /* always destroy the old matrix and create a new memory; 129 hope this does not churn the memory too much. The problem 130 is I do not know if it is possible to put the matrix back to 131 its initial state so that we can directly copy the values 132 the second time through. */ 133 PetscStackCallStandard(HYPRE_IJMatrixDestroy,(jac->ij)); 134 jac->ij = 0; 135 } 136 137 if (!jac->ij) { /* create the matrix the first time through */ 138 ierr = MatHYPRE_IJMatrixCreate(pc->pmat,&jac->ij);CHKERRQ(ierr); 139 } 140 if (!jac->b) { /* create the vectors the first time through */ 141 Vec x,b; 142 ierr = MatCreateVecs(pc->pmat,&x,&b);CHKERRQ(ierr); 143 ierr = VecHYPRE_IJVectorCreate(x,&jac->x);CHKERRQ(ierr); 144 ierr = VecHYPRE_IJVectorCreate(b,&jac->b);CHKERRQ(ierr); 145 ierr = VecDestroy(&x);CHKERRQ(ierr); 146 ierr = VecDestroy(&b);CHKERRQ(ierr); 147 } 148 149 /* special case for BoomerAMG */ 150 if (jac->setup == HYPRE_BoomerAMGSetup) { 151 ierr = MatGetBlockSize(pc->pmat,&bs);CHKERRQ(ierr); 152 if (bs > 1) PetscStackCallStandard(HYPRE_BoomerAMGSetNumFunctions,(jac->hsolver,bs)); 153 } 154 155 /* special case for AMS */ 156 if (jac->setup == HYPRE_AMSSetup) { 157 if (!jac->coords[0] && !jac->constants[0]) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE AMS preconditioner needs either coordinate vectors via PCSetCoordinates() or edge constant vectors via PCHYPRESetEdgeConstantVectors()"); 158 if (!jac->G) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE AMS preconditioner needs discrete gradient operator via PCHYPRESetDiscreteGradient"); 159 } 160 /* special case for ADS */ 161 if (jac->setup == HYPRE_ADSSetup) { 162 if (!jac->coords[0]) { 163 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs coordinate vectors via PCSetCoordinates()"); 164 } else if (!jac->coords[1] || !jac->coords[2]) { 165 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner has been designed for three dimensional problems! For two dimensional problems, use HYPRE AMS instead"); 166 } 167 if (!jac->G) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs discrete gradient operator via PCHYPRESetDiscreteGradient"); 168 if (!jac->C) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs discrete curl operator via PCHYPRESetDiscreteGradient"); 169 } 170 ierr = MatHYPRE_IJMatrixCopy(pc->pmat,jac->ij);CHKERRQ(ierr); 171 PetscStackCallStandard(HYPRE_IJMatrixGetObject,(jac->ij,(void**)&hmat)); 172 PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->b,(void**)&bv)); 173 PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->x,(void**)&xv)); 174 PetscStackCall("HYPRE_SetupXXX",ierr = (*jac->setup)(jac->hsolver,hmat,bv,xv);CHKERRQ(ierr);); 175 PetscFunctionReturn(0); 176 } 177 178 /* 179 Replaces the address where the HYPRE vector points to its data with the address of 180 PETSc's data. Saves the old address so it can be reset when we are finished with it. 181 Allows use to get the data into a HYPRE vector without the cost of memcopies 182 */ 183 #define HYPREReplacePointer(b,newvalue,savedvalue) { \ 184 hypre_ParVector *par_vector = (hypre_ParVector*)hypre_IJVectorObject(((hypre_IJVector*)b)); \ 185 hypre_Vector *local_vector = hypre_ParVectorLocalVector(par_vector); \ 186 savedvalue = local_vector->data; \ 187 local_vector->data = newvalue; \ 188 } 189 190 #undef __FUNCT__ 191 #define __FUNCT__ "PCApply_HYPRE" 192 static PetscErrorCode PCApply_HYPRE(PC pc,Vec b,Vec x) 193 { 194 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 195 PetscErrorCode ierr; 196 HYPRE_ParCSRMatrix hmat; 197 PetscScalar *xv; 198 const PetscScalar *bv,*sbv; 199 HYPRE_ParVector jbv,jxv; 200 PetscScalar *sxv; 201 PetscInt hierr; 202 203 PetscFunctionBegin; 204 ierr = PetscCitationsRegister(hypreCitation,&cite);CHKERRQ(ierr); 205 if (!jac->applyrichardson) {ierr = VecSet(x,0.0);CHKERRQ(ierr);} 206 ierr = VecGetArrayRead(b,&bv);CHKERRQ(ierr); 207 ierr = VecGetArray(x,&xv);CHKERRQ(ierr); 208 HYPREReplacePointer(jac->b,(PetscScalar*)bv,sbv); 209 HYPREReplacePointer(jac->x,xv,sxv); 210 211 PetscStackCallStandard(HYPRE_IJMatrixGetObject,(jac->ij,(void**)&hmat)); 212 PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->b,(void**)&jbv)); 213 PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->x,(void**)&jxv)); 214 PetscStackCall("Hypre solve",hierr = (*jac->solve)(jac->hsolver,hmat,jbv,jxv); 215 if (hierr && hierr != HYPRE_ERROR_CONV) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",hierr); 216 if (hierr) hypre__global_error = 0;); 217 218 HYPREReplacePointer(jac->b,(PetscScalar*)sbv,bv); 219 HYPREReplacePointer(jac->x,sxv,xv); 220 ierr = VecRestoreArray(x,&xv);CHKERRQ(ierr); 221 ierr = VecRestoreArrayRead(b,&bv);CHKERRQ(ierr); 222 PetscFunctionReturn(0); 223 } 224 225 #undef __FUNCT__ 226 #define __FUNCT__ "PCDestroy_HYPRE" 227 static PetscErrorCode PCDestroy_HYPRE(PC pc) 228 { 229 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 230 PetscErrorCode ierr; 231 232 PetscFunctionBegin; 233 if (jac->ij) PetscStackCallStandard(HYPRE_IJMatrixDestroy,(jac->ij)); 234 if (jac->b) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->b)); 235 if (jac->x) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->x)); 236 if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[0])); 237 if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[1])); 238 if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[2])); 239 if (jac->constants[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[0])); 240 if (jac->constants[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[1])); 241 if (jac->constants[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[2])); 242 if (jac->G) PetscStackCallStandard(HYPRE_IJMatrixDestroy,(jac->G)); 243 if (jac->C) PetscStackCallStandard(HYPRE_IJMatrixDestroy,(jac->C)); 244 if (jac->alpha_Poisson) PetscStackCallStandard(HYPRE_IJMatrixDestroy,(jac->alpha_Poisson)); 245 if (jac->beta_Poisson) PetscStackCallStandard(HYPRE_IJMatrixDestroy,(jac->beta_Poisson)); 246 if (jac->destroy) PetscStackCall("HYPRE_DestroyXXX",ierr = (*jac->destroy)(jac->hsolver);CHKERRQ(ierr);); 247 ierr = PetscFree(jac->hypre_type);CHKERRQ(ierr); 248 if (jac->comm_hypre != MPI_COMM_NULL) { ierr = MPI_Comm_free(&(jac->comm_hypre));CHKERRQ(ierr);} 249 ierr = PetscFree(pc->data);CHKERRQ(ierr); 250 251 ierr = PetscObjectChangeTypeName((PetscObject)pc,0);CHKERRQ(ierr); 252 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetType_C",NULL);CHKERRQ(ierr); 253 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPREGetType_C",NULL);CHKERRQ(ierr); 254 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetCoordinates_C",NULL);CHKERRQ(ierr); 255 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteGradient_C",NULL);CHKERRQ(ierr); 256 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteCurl_C",NULL);CHKERRQ(ierr); 257 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetConstantEdgeVectors_C",NULL);CHKERRQ(ierr); 258 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetAlphaPoissonMatrix_C",NULL);CHKERRQ(ierr); 259 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetBetaPoissonMatrix_C",NULL);CHKERRQ(ierr); 260 PetscFunctionReturn(0); 261 } 262 263 /* --------------------------------------------------------------------------------------------*/ 264 #undef __FUNCT__ 265 #define __FUNCT__ "PCSetFromOptions_HYPRE_Pilut" 266 static PetscErrorCode PCSetFromOptions_HYPRE_Pilut(PetscOptions *PetscOptionsObject,PC pc) 267 { 268 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 269 PetscErrorCode ierr; 270 PetscBool flag; 271 272 PetscFunctionBegin; 273 ierr = PetscOptionsHead(PetscOptionsObject,"HYPRE Pilut Options");CHKERRQ(ierr); 274 ierr = PetscOptionsInt("-pc_hypre_pilut_maxiter","Number of iterations","None",jac->maxiter,&jac->maxiter,&flag);CHKERRQ(ierr); 275 if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetMaxIter,(jac->hsolver,jac->maxiter)); 276 ierr = PetscOptionsReal("-pc_hypre_pilut_tol","Drop tolerance","None",jac->tol,&jac->tol,&flag);CHKERRQ(ierr); 277 if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetDropTolerance,(jac->hsolver,jac->tol)); 278 ierr = PetscOptionsInt("-pc_hypre_pilut_factorrowsize","FactorRowSize","None",jac->factorrowsize,&jac->factorrowsize,&flag);CHKERRQ(ierr); 279 if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetFactorRowSize,(jac->hsolver,jac->factorrowsize)); 280 ierr = PetscOptionsTail();CHKERRQ(ierr); 281 PetscFunctionReturn(0); 282 } 283 284 #undef __FUNCT__ 285 #define __FUNCT__ "PCView_HYPRE_Pilut" 286 static PetscErrorCode PCView_HYPRE_Pilut(PC pc,PetscViewer viewer) 287 { 288 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 289 PetscErrorCode ierr; 290 PetscBool iascii; 291 292 PetscFunctionBegin; 293 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 294 if (iascii) { 295 ierr = PetscViewerASCIIPrintf(viewer," HYPRE Pilut preconditioning\n");CHKERRQ(ierr); 296 if (jac->maxiter != PETSC_DEFAULT) { 297 ierr = PetscViewerASCIIPrintf(viewer," HYPRE Pilut: maximum number of iterations %d\n",jac->maxiter);CHKERRQ(ierr); 298 } else { 299 ierr = PetscViewerASCIIPrintf(viewer," HYPRE Pilut: default maximum number of iterations \n");CHKERRQ(ierr); 300 } 301 if (jac->tol != PETSC_DEFAULT) { 302 ierr = PetscViewerASCIIPrintf(viewer," HYPRE Pilut: drop tolerance %g\n",(double)jac->tol);CHKERRQ(ierr); 303 } else { 304 ierr = PetscViewerASCIIPrintf(viewer," HYPRE Pilut: default drop tolerance \n");CHKERRQ(ierr); 305 } 306 if (jac->factorrowsize != PETSC_DEFAULT) { 307 ierr = PetscViewerASCIIPrintf(viewer," HYPRE Pilut: factor row size %d\n",jac->factorrowsize);CHKERRQ(ierr); 308 } else { 309 ierr = PetscViewerASCIIPrintf(viewer," HYPRE Pilut: default factor row size \n");CHKERRQ(ierr); 310 } 311 } 312 PetscFunctionReturn(0); 313 } 314 315 /* --------------------------------------------------------------------------------------------*/ 316 317 #undef __FUNCT__ 318 #define __FUNCT__ "PCApplyTranspose_HYPRE_BoomerAMG" 319 static PetscErrorCode PCApplyTranspose_HYPRE_BoomerAMG(PC pc,Vec b,Vec x) 320 { 321 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 322 PetscErrorCode ierr; 323 HYPRE_ParCSRMatrix hmat; 324 PetscScalar *xv; 325 const PetscScalar *bv; 326 HYPRE_ParVector jbv,jxv; 327 PetscScalar *sbv,*sxv; 328 PetscInt hierr; 329 330 PetscFunctionBegin; 331 ierr = PetscCitationsRegister(hypreCitation,&cite);CHKERRQ(ierr); 332 ierr = VecSet(x,0.0);CHKERRQ(ierr); 333 ierr = VecGetArrayRead(b,&bv);CHKERRQ(ierr); 334 ierr = VecGetArray(x,&xv);CHKERRQ(ierr); 335 HYPREReplacePointer(jac->b,(PetscScalar*)bv,sbv); 336 HYPREReplacePointer(jac->x,xv,sxv); 337 338 PetscStackCallStandard(HYPRE_IJMatrixGetObject,(jac->ij,(void**)&hmat)); 339 PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->b,(void**)&jbv)); 340 PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->x,(void**)&jxv)); 341 342 hierr = HYPRE_BoomerAMGSolveT(jac->hsolver,hmat,jbv,jxv); 343 /* error code of 1 in BoomerAMG merely means convergence not achieved */ 344 if (hierr && (hierr != 1)) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",hierr); 345 if (hierr) hypre__global_error = 0; 346 347 HYPREReplacePointer(jac->b,sbv,bv); 348 HYPREReplacePointer(jac->x,sxv,xv); 349 ierr = VecRestoreArray(x,&xv);CHKERRQ(ierr); 350 ierr = VecRestoreArrayRead(b,&bv);CHKERRQ(ierr); 351 PetscFunctionReturn(0); 352 } 353 354 /* static array length */ 355 #define ALEN(a) (sizeof(a)/sizeof((a)[0])) 356 357 static const char *HYPREBoomerAMGCycleType[] = {"","V","W"}; 358 static const char *HYPREBoomerAMGCoarsenType[] = {"CLJP","Ruge-Stueben","","modifiedRuge-Stueben","","","Falgout", "", "PMIS", "", "HMIS"}; 359 static const char *HYPREBoomerAMGMeasureType[] = {"local","global"}; 360 /* The following corresponds to HYPRE_BoomerAMGSetRelaxType which has many missing numbers in the enum */ 361 static const char *HYPREBoomerAMGRelaxType[] = {"Jacobi","sequential-Gauss-Seidel","seqboundary-Gauss-Seidel","SOR/Jacobi","backward-SOR/Jacobi", 362 "" /* [5] hybrid chaotic Gauss-Seidel (works only with OpenMP) */,"symmetric-SOR/Jacobi", 363 "" /* 7 */,"l1scaled-SOR/Jacobi","Gaussian-elimination", 364 "" /* 10 */, "" /* 11 */, "" /* 12 */, "" /* 13 */, "" /* 14 */, 365 "CG" /* non-stationary */,"Chebyshev","FCF-Jacobi","l1scaled-Jacobi"}; 366 static const char *HYPREBoomerAMGInterpType[] = {"classical", "", "", "direct", "multipass", "multipass-wts", "ext+i", 367 "ext+i-cc", "standard", "standard-wts", "", "", "FF", "FF1"}; 368 #undef __FUNCT__ 369 #define __FUNCT__ "PCSetFromOptions_HYPRE_BoomerAMG" 370 static PetscErrorCode PCSetFromOptions_HYPRE_BoomerAMG(PetscOptions *PetscOptionsObject,PC pc) 371 { 372 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 373 PetscErrorCode ierr; 374 PetscInt n,indx,level; 375 PetscBool flg, tmp_truth; 376 double tmpdbl, twodbl[2]; 377 378 PetscFunctionBegin; 379 ierr = PetscOptionsHead(PetscOptionsObject,"HYPRE BoomerAMG Options");CHKERRQ(ierr); 380 ierr = PetscOptionsEList("-pc_hypre_boomeramg_cycle_type","Cycle type","None",HYPREBoomerAMGCycleType+1,2,HYPREBoomerAMGCycleType[jac->cycletype],&indx,&flg);CHKERRQ(ierr); 381 if (flg) { 382 jac->cycletype = indx+1; 383 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleType,(jac->hsolver,jac->cycletype)); 384 } 385 ierr = PetscOptionsInt("-pc_hypre_boomeramg_max_levels","Number of levels (of grids) allowed","None",jac->maxlevels,&jac->maxlevels,&flg);CHKERRQ(ierr); 386 if (flg) { 387 if (jac->maxlevels < 2) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of levels %d must be at least two",jac->maxlevels); 388 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxLevels,(jac->hsolver,jac->maxlevels)); 389 } 390 ierr = PetscOptionsInt("-pc_hypre_boomeramg_max_iter","Maximum iterations used PER hypre call","None",jac->maxiter,&jac->maxiter,&flg);CHKERRQ(ierr); 391 if (flg) { 392 if (jac->maxiter < 1) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of iterations %d must be at least one",jac->maxiter); 393 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,jac->maxiter)); 394 } 395 ierr = PetscOptionsScalar("-pc_hypre_boomeramg_tol","Convergence tolerance PER hypre call (0.0 = use a fixed number of iterations)","None",jac->tol,&jac->tol,&flg);CHKERRQ(ierr); 396 if (flg) { 397 if (jac->tol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Tolerance %g must be greater than or equal to zero",(double)jac->tol); 398 PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,jac->tol)); 399 } 400 401 ierr = PetscOptionsScalar("-pc_hypre_boomeramg_truncfactor","Truncation factor for interpolation (0=no truncation)","None",jac->truncfactor,&jac->truncfactor,&flg);CHKERRQ(ierr); 402 if (flg) { 403 if (jac->truncfactor < 0.0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Truncation factor %g must be great than or equal zero",(double)jac->truncfactor); 404 PetscStackCallStandard(HYPRE_BoomerAMGSetTruncFactor,(jac->hsolver,jac->truncfactor)); 405 } 406 407 ierr = PetscOptionsInt("-pc_hypre_boomeramg_P_max","Max elements per row for interpolation operator (0=unlimited)","None",jac->pmax,&jac->pmax,&flg);CHKERRQ(ierr); 408 if (flg) { 409 if (jac->pmax < 0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"P_max %g must be greater than or equal to zero",(double)jac->pmax); 410 PetscStackCallStandard(HYPRE_BoomerAMGSetPMaxElmts,(jac->hsolver,jac->pmax)); 411 } 412 413 ierr = PetscOptionsInt("-pc_hypre_boomeramg_agg_nl","Number of levels of aggressive coarsening","None",jac->agg_nl,&jac->agg_nl,&flg);CHKERRQ(ierr); 414 if (flg) { 415 if (jac->agg_nl < 0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of levels %g must be greater than or equal to zero",(double)jac->agg_nl); 416 417 PetscStackCallStandard(HYPRE_BoomerAMGSetAggNumLevels,(jac->hsolver,jac->agg_nl)); 418 } 419 420 421 ierr = PetscOptionsInt("-pc_hypre_boomeramg_agg_num_paths","Number of paths for aggressive coarsening","None",jac->agg_num_paths,&jac->agg_num_paths,&flg);CHKERRQ(ierr); 422 if (flg) { 423 if (jac->agg_num_paths < 1) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of paths %g must be greater than or equal to 1",(double)jac->agg_num_paths); 424 425 PetscStackCallStandard(HYPRE_BoomerAMGSetNumPaths,(jac->hsolver,jac->agg_num_paths)); 426 } 427 428 429 ierr = PetscOptionsScalar("-pc_hypre_boomeramg_strong_threshold","Threshold for being strongly connected","None",jac->strongthreshold,&jac->strongthreshold,&flg);CHKERRQ(ierr); 430 if (flg) { 431 if (jac->strongthreshold < 0.0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Strong threshold %g must be great than or equal zero",(double)jac->strongthreshold); 432 PetscStackCallStandard(HYPRE_BoomerAMGSetStrongThreshold,(jac->hsolver,jac->strongthreshold)); 433 } 434 ierr = PetscOptionsScalar("-pc_hypre_boomeramg_max_row_sum","Maximum row sum","None",jac->maxrowsum,&jac->maxrowsum,&flg);CHKERRQ(ierr); 435 if (flg) { 436 if (jac->maxrowsum < 0.0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Maximum row sum %g must be greater than zero",(double)jac->maxrowsum); 437 if (jac->maxrowsum > 1.0) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Maximum row sum %g must be less than or equal one",(double)jac->maxrowsum); 438 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxRowSum,(jac->hsolver,jac->maxrowsum)); 439 } 440 441 /* Grid sweeps */ 442 ierr = PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_all","Number of sweeps for the up and down grid levels","None",jac->gridsweeps[0],&indx,&flg);CHKERRQ(ierr); 443 if (flg) { 444 PetscStackCallStandard(HYPRE_BoomerAMGSetNumSweeps,(jac->hsolver,indx)); 445 /* modify the jac structure so we can view the updated options with PC_View */ 446 jac->gridsweeps[0] = indx; 447 jac->gridsweeps[1] = indx; 448 /*defaults coarse to 1 */ 449 jac->gridsweeps[2] = 1; 450 } 451 452 ierr = PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_down","Number of sweeps for the down cycles","None",jac->gridsweeps[0], &indx,&flg);CHKERRQ(ierr); 453 if (flg) { 454 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,(jac->hsolver,indx, 1)); 455 jac->gridsweeps[0] = indx; 456 } 457 ierr = PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_up","Number of sweeps for the up cycles","None",jac->gridsweeps[1],&indx,&flg);CHKERRQ(ierr); 458 if (flg) { 459 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,(jac->hsolver,indx, 2)); 460 jac->gridsweeps[1] = indx; 461 } 462 ierr = PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_coarse","Number of sweeps for the coarse level","None",jac->gridsweeps[2],&indx,&flg);CHKERRQ(ierr); 463 if (flg) { 464 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,(jac->hsolver,indx, 3)); 465 jac->gridsweeps[2] = indx; 466 } 467 468 /* Relax type */ 469 ierr = PetscOptionsEList("-pc_hypre_boomeramg_relax_type_all","Relax type for the up and down cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg);CHKERRQ(ierr); 470 if (flg) { 471 jac->relaxtype[0] = jac->relaxtype[1] = indx; 472 PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxType,(jac->hsolver, indx)); 473 /* by default, coarse type set to 9 */ 474 jac->relaxtype[2] = 9; 475 476 } 477 ierr = PetscOptionsEList("-pc_hypre_boomeramg_relax_type_down","Relax type for the down cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg);CHKERRQ(ierr); 478 if (flg) { 479 jac->relaxtype[0] = indx; 480 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,(jac->hsolver, indx, 1)); 481 } 482 ierr = PetscOptionsEList("-pc_hypre_boomeramg_relax_type_up","Relax type for the up cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg);CHKERRQ(ierr); 483 if (flg) { 484 jac->relaxtype[1] = indx; 485 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,(jac->hsolver, indx, 2)); 486 } 487 ierr = PetscOptionsEList("-pc_hypre_boomeramg_relax_type_coarse","Relax type on coarse grid","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[9],&indx,&flg);CHKERRQ(ierr); 488 if (flg) { 489 jac->relaxtype[2] = indx; 490 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,(jac->hsolver, indx, 3)); 491 } 492 493 /* Relaxation Weight */ 494 ierr = PetscOptionsReal("-pc_hypre_boomeramg_relax_weight_all","Relaxation weight for all levels (0 = hypre estimates, -k = determined with k CG steps)","None",jac->relaxweight, &tmpdbl,&flg);CHKERRQ(ierr); 495 if (flg) { 496 PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxWt,(jac->hsolver,tmpdbl)); 497 jac->relaxweight = tmpdbl; 498 } 499 500 n = 2; 501 twodbl[0] = twodbl[1] = 1.0; 502 ierr = PetscOptionsRealArray("-pc_hypre_boomeramg_relax_weight_level","Set the relaxation weight for a particular level (weight,level)","None",twodbl, &n, &flg);CHKERRQ(ierr); 503 if (flg) { 504 if (n == 2) { 505 indx = (int)PetscAbsReal(twodbl[1]); 506 PetscStackCallStandard(HYPRE_BoomerAMGSetLevelRelaxWt,(jac->hsolver,twodbl[0],indx)); 507 } else SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Relax weight level: you must provide 2 values separated by a comma (and no space), you provided %d",n); 508 } 509 510 /* Outer relaxation Weight */ 511 ierr = PetscOptionsReal("-pc_hypre_boomeramg_outer_relax_weight_all","Outer relaxation weight for all levels (-k = determined with k CG steps)","None",jac->outerrelaxweight, &tmpdbl,&flg);CHKERRQ(ierr); 512 if (flg) { 513 PetscStackCallStandard(HYPRE_BoomerAMGSetOuterWt,(jac->hsolver, tmpdbl)); 514 jac->outerrelaxweight = tmpdbl; 515 } 516 517 n = 2; 518 twodbl[0] = twodbl[1] = 1.0; 519 ierr = PetscOptionsRealArray("-pc_hypre_boomeramg_outer_relax_weight_level","Set the outer relaxation weight for a particular level (weight,level)","None",twodbl, &n, &flg);CHKERRQ(ierr); 520 if (flg) { 521 if (n == 2) { 522 indx = (int)PetscAbsReal(twodbl[1]); 523 PetscStackCallStandard(HYPRE_BoomerAMGSetLevelOuterWt,(jac->hsolver, twodbl[0], indx)); 524 } else SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Relax weight outer level: You must provide 2 values separated by a comma (and no space), you provided %d",n); 525 } 526 527 /* the Relax Order */ 528 ierr = PetscOptionsBool("-pc_hypre_boomeramg_no_CF", "Do not use CF-relaxation", "None", PETSC_FALSE, &tmp_truth, &flg);CHKERRQ(ierr); 529 530 if (flg && tmp_truth) { 531 jac->relaxorder = 0; 532 PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxOrder,(jac->hsolver, jac->relaxorder)); 533 } 534 ierr = PetscOptionsEList("-pc_hypre_boomeramg_measure_type","Measure type","None",HYPREBoomerAMGMeasureType,ALEN(HYPREBoomerAMGMeasureType),HYPREBoomerAMGMeasureType[0],&indx,&flg);CHKERRQ(ierr); 535 if (flg) { 536 jac->measuretype = indx; 537 PetscStackCallStandard(HYPRE_BoomerAMGSetMeasureType,(jac->hsolver,jac->measuretype)); 538 } 539 /* update list length 3/07 */ 540 ierr = PetscOptionsEList("-pc_hypre_boomeramg_coarsen_type","Coarsen type","None",HYPREBoomerAMGCoarsenType,ALEN(HYPREBoomerAMGCoarsenType),HYPREBoomerAMGCoarsenType[6],&indx,&flg);CHKERRQ(ierr); 541 if (flg) { 542 jac->coarsentype = indx; 543 PetscStackCallStandard(HYPRE_BoomerAMGSetCoarsenType,(jac->hsolver,jac->coarsentype)); 544 } 545 546 /* new 3/07 */ 547 ierr = PetscOptionsEList("-pc_hypre_boomeramg_interp_type","Interpolation type","None",HYPREBoomerAMGInterpType,ALEN(HYPREBoomerAMGInterpType),HYPREBoomerAMGInterpType[0],&indx,&flg);CHKERRQ(ierr); 548 if (flg) { 549 jac->interptype = indx; 550 PetscStackCallStandard(HYPRE_BoomerAMGSetInterpType,(jac->hsolver,jac->interptype)); 551 } 552 553 ierr = PetscOptionsName("-pc_hypre_boomeramg_print_statistics","Print statistics","None",&flg);CHKERRQ(ierr); 554 if (flg) { 555 level = 3; 556 ierr = PetscOptionsInt("-pc_hypre_boomeramg_print_statistics","Print statistics","None",level,&level,NULL);CHKERRQ(ierr); 557 558 jac->printstatistics = PETSC_TRUE; 559 PetscStackCallStandard(HYPRE_BoomerAMGSetPrintLevel,(jac->hsolver,level)); 560 } 561 562 ierr = PetscOptionsName("-pc_hypre_boomeramg_print_debug","Print debug information","None",&flg);CHKERRQ(ierr); 563 if (flg) { 564 level = 3; 565 ierr = PetscOptionsInt("-pc_hypre_boomeramg_print_debug","Print debug information","None",level,&level,NULL);CHKERRQ(ierr); 566 567 jac->printstatistics = PETSC_TRUE; 568 PetscStackCallStandard(HYPRE_BoomerAMGSetDebugFlag,(jac->hsolver,level)); 569 } 570 571 ierr = PetscOptionsBool("-pc_hypre_boomeramg_nodal_coarsen", "HYPRE_BoomerAMGSetNodal()", "None", jac->nodal_coarsen ? PETSC_TRUE : PETSC_FALSE, &tmp_truth, &flg);CHKERRQ(ierr); 572 if (flg && tmp_truth) { 573 jac->nodal_coarsen = 1; 574 PetscStackCallStandard(HYPRE_BoomerAMGSetNodal,(jac->hsolver,1)); 575 } 576 577 ierr = PetscOptionsBool("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None", PETSC_FALSE, &tmp_truth, &flg);CHKERRQ(ierr); 578 if (flg && tmp_truth) { 579 PetscInt tmp_int; 580 ierr = PetscOptionsInt("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None",jac->nodal_relax_levels,&tmp_int,&flg);CHKERRQ(ierr); 581 if (flg) jac->nodal_relax_levels = tmp_int; 582 PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothType,(jac->hsolver,6)); 583 PetscStackCallStandard(HYPRE_BoomerAMGSetDomainType,(jac->hsolver,1)); 584 PetscStackCallStandard(HYPRE_BoomerAMGSetOverlap,(jac->hsolver,0)); 585 PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothNumLevels,(jac->hsolver,jac->nodal_relax_levels)); 586 } 587 588 ierr = PetscOptionsTail();CHKERRQ(ierr); 589 PetscFunctionReturn(0); 590 } 591 592 #undef __FUNCT__ 593 #define __FUNCT__ "PCApplyRichardson_HYPRE_BoomerAMG" 594 static PetscErrorCode PCApplyRichardson_HYPRE_BoomerAMG(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscBool guesszero,PetscInt *outits,PCRichardsonConvergedReason *reason) 595 { 596 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 597 PetscErrorCode ierr; 598 PetscInt oits; 599 600 PetscFunctionBegin; 601 ierr = PetscCitationsRegister(hypreCitation,&cite);CHKERRQ(ierr); 602 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,its*jac->maxiter)); 603 PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,rtol)); 604 jac->applyrichardson = PETSC_TRUE; 605 ierr = PCApply_HYPRE(pc,b,y);CHKERRQ(ierr); 606 jac->applyrichardson = PETSC_FALSE; 607 PetscStackCallStandard(HYPRE_BoomerAMGGetNumIterations,(jac->hsolver,(HYPRE_Int *)&oits)); 608 *outits = oits; 609 if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS; 610 else *reason = PCRICHARDSON_CONVERGED_RTOL; 611 PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,jac->tol)); 612 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,jac->maxiter)); 613 PetscFunctionReturn(0); 614 } 615 616 617 #undef __FUNCT__ 618 #define __FUNCT__ "PCView_HYPRE_BoomerAMG" 619 static PetscErrorCode PCView_HYPRE_BoomerAMG(PC pc,PetscViewer viewer) 620 { 621 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 622 PetscErrorCode ierr; 623 PetscBool iascii; 624 625 PetscFunctionBegin; 626 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 627 if (iascii) { 628 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG preconditioning\n");CHKERRQ(ierr); 629 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Cycle type %s\n",HYPREBoomerAMGCycleType[jac->cycletype]);CHKERRQ(ierr); 630 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Maximum number of levels %d\n",jac->maxlevels);CHKERRQ(ierr); 631 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Maximum number of iterations PER hypre call %d\n",jac->maxiter);CHKERRQ(ierr); 632 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Convergence tolerance PER hypre call %g\n",(double)jac->tol);CHKERRQ(ierr); 633 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Threshold for strong coupling %g\n",(double)jac->strongthreshold);CHKERRQ(ierr); 634 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Interpolation truncation factor %g\n",(double)jac->truncfactor);CHKERRQ(ierr); 635 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Interpolation: max elements per row %d\n",jac->pmax);CHKERRQ(ierr); 636 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Number of levels of aggressive coarsening %d\n",jac->agg_nl);CHKERRQ(ierr); 637 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Number of paths for aggressive coarsening %d\n",jac->agg_num_paths);CHKERRQ(ierr); 638 639 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Maximum row sums %g\n",(double)jac->maxrowsum);CHKERRQ(ierr); 640 641 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Sweeps down %d\n",jac->gridsweeps[0]);CHKERRQ(ierr); 642 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Sweeps up %d\n",jac->gridsweeps[1]);CHKERRQ(ierr); 643 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Sweeps on coarse %d\n",jac->gridsweeps[2]);CHKERRQ(ierr); 644 645 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Relax down %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[0]]);CHKERRQ(ierr); 646 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Relax up %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[1]]);CHKERRQ(ierr); 647 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Relax on coarse %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[2]]);CHKERRQ(ierr); 648 649 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Relax weight (all) %g\n",(double)jac->relaxweight);CHKERRQ(ierr); 650 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Outer relax weight (all) %g\n",(double)jac->outerrelaxweight);CHKERRQ(ierr); 651 652 if (jac->relaxorder) { 653 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Using CF-relaxation\n");CHKERRQ(ierr); 654 } else { 655 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Not using CF-relaxation\n");CHKERRQ(ierr); 656 } 657 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Measure type %s\n",HYPREBoomerAMGMeasureType[jac->measuretype]);CHKERRQ(ierr); 658 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Coarsen type %s\n",HYPREBoomerAMGCoarsenType[jac->coarsentype]);CHKERRQ(ierr); 659 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Interpolation type %s\n",HYPREBoomerAMGInterpType[jac->interptype]);CHKERRQ(ierr); 660 if (jac->nodal_coarsen) { 661 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Using nodal coarsening (with HYPRE_BOOMERAMGSetNodal())\n");CHKERRQ(ierr); 662 } 663 if (jac->nodal_relax) { 664 ierr = PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Using nodal relaxation via Schwarz smoothing on levels %d\n",jac->nodal_relax_levels);CHKERRQ(ierr); 665 } 666 } 667 PetscFunctionReturn(0); 668 } 669 670 /* --------------------------------------------------------------------------------------------*/ 671 #undef __FUNCT__ 672 #define __FUNCT__ "PCSetFromOptions_HYPRE_ParaSails" 673 static PetscErrorCode PCSetFromOptions_HYPRE_ParaSails(PetscOptions *PetscOptionsObject,PC pc) 674 { 675 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 676 PetscErrorCode ierr; 677 PetscInt indx; 678 PetscBool flag; 679 const char *symtlist[] = {"nonsymmetric","SPD","nonsymmetric,SPD"}; 680 681 PetscFunctionBegin; 682 ierr = PetscOptionsHead(PetscOptionsObject,"HYPRE ParaSails Options");CHKERRQ(ierr); 683 ierr = PetscOptionsInt("-pc_hypre_parasails_nlevels","Number of number of levels","None",jac->nlevels,&jac->nlevels,0);CHKERRQ(ierr); 684 ierr = PetscOptionsReal("-pc_hypre_parasails_thresh","Threshold","None",jac->threshhold,&jac->threshhold,&flag);CHKERRQ(ierr); 685 if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetParams,(jac->hsolver,jac->threshhold,jac->nlevels)); 686 687 ierr = PetscOptionsReal("-pc_hypre_parasails_filter","filter","None",jac->filter,&jac->filter,&flag);CHKERRQ(ierr); 688 if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetFilter,(jac->hsolver,jac->filter)); 689 690 ierr = PetscOptionsReal("-pc_hypre_parasails_loadbal","Load balance","None",jac->loadbal,&jac->loadbal,&flag);CHKERRQ(ierr); 691 if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetLoadbal,(jac->hsolver,jac->loadbal)); 692 693 ierr = PetscOptionsBool("-pc_hypre_parasails_logging","Print info to screen","None",(PetscBool)jac->logging,(PetscBool*)&jac->logging,&flag);CHKERRQ(ierr); 694 if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetLogging,(jac->hsolver,jac->logging)); 695 696 ierr = PetscOptionsBool("-pc_hypre_parasails_reuse","Reuse nonzero pattern in preconditioner","None",(PetscBool)jac->ruse,(PetscBool*)&jac->ruse,&flag);CHKERRQ(ierr); 697 if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetReuse,(jac->hsolver,jac->ruse)); 698 699 ierr = PetscOptionsEList("-pc_hypre_parasails_sym","Symmetry of matrix and preconditioner","None",symtlist,ALEN(symtlist),symtlist[0],&indx,&flag);CHKERRQ(ierr); 700 if (flag) { 701 jac->symt = indx; 702 PetscStackCallStandard(HYPRE_ParaSailsSetSym,(jac->hsolver,jac->symt)); 703 } 704 705 ierr = PetscOptionsTail();CHKERRQ(ierr); 706 PetscFunctionReturn(0); 707 } 708 709 #undef __FUNCT__ 710 #define __FUNCT__ "PCView_HYPRE_ParaSails" 711 static PetscErrorCode PCView_HYPRE_ParaSails(PC pc,PetscViewer viewer) 712 { 713 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 714 PetscErrorCode ierr; 715 PetscBool iascii; 716 const char *symt = 0;; 717 718 PetscFunctionBegin; 719 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 720 if (iascii) { 721 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ParaSails preconditioning\n");CHKERRQ(ierr); 722 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: nlevels %d\n",jac->nlevels);CHKERRQ(ierr); 723 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: threshold %g\n",(double)jac->threshhold);CHKERRQ(ierr); 724 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: filter %g\n",(double)jac->filter);CHKERRQ(ierr); 725 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: load balance %g\n",(double)jac->loadbal);CHKERRQ(ierr); 726 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: reuse nonzero structure %s\n",PetscBools[jac->ruse]);CHKERRQ(ierr); 727 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: print info to screen %s\n",PetscBools[jac->logging]);CHKERRQ(ierr); 728 if (!jac->symt) symt = "nonsymmetric matrix and preconditioner"; 729 else if (jac->symt == 1) symt = "SPD matrix and preconditioner"; 730 else if (jac->symt == 2) symt = "nonsymmetric matrix but SPD preconditioner"; 731 else SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG,"Unknown HYPRE ParaSails symmetric option %d",jac->symt); 732 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: %s\n",symt);CHKERRQ(ierr); 733 } 734 PetscFunctionReturn(0); 735 } 736 /* --------------------------------------------------------------------------------------------*/ 737 #undef __FUNCT__ 738 #define __FUNCT__ "PCSetFromOptions_HYPRE_AMS" 739 static PetscErrorCode PCSetFromOptions_HYPRE_AMS(PetscOptions *PetscOptionsObject,PC pc) 740 { 741 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 742 PetscErrorCode ierr; 743 PetscInt n; 744 PetscBool flag,flag2,flag3,flag4; 745 746 PetscFunctionBegin; 747 ierr = PetscOptionsHead(PetscOptionsObject,"HYPRE AMS Options");CHKERRQ(ierr); 748 ierr = PetscOptionsInt("-pc_hypre_ams_print_level","Debugging output level for AMS","None",jac->as_print,&jac->as_print,&flag);CHKERRQ(ierr); 749 if (flag) PetscStackCallStandard(HYPRE_AMSSetPrintLevel,(jac->hsolver,jac->as_print)); 750 ierr = PetscOptionsInt("-pc_hypre_ams_max_iter","Maximum number of AMS multigrid iterations within PCApply","None",jac->as_max_iter,&jac->as_max_iter,&flag);CHKERRQ(ierr); 751 if (flag) PetscStackCallStandard(HYPRE_AMSSetMaxIter,(jac->hsolver,jac->as_max_iter)); 752 ierr = PetscOptionsInt("-pc_hypre_ams_cycle_type","Cycle type for AMS multigrid","None",jac->ams_cycle_type,&jac->ams_cycle_type,&flag);CHKERRQ(ierr); 753 if (flag) PetscStackCallStandard(HYPRE_AMSSetCycleType,(jac->hsolver,jac->ams_cycle_type)); 754 ierr = PetscOptionsReal("-pc_hypre_ams_tol","Error tolerance for AMS multigrid","None",jac->as_tol,&jac->as_tol,&flag);CHKERRQ(ierr); 755 if (flag) PetscStackCallStandard(HYPRE_AMSSetTol,(jac->hsolver,jac->as_tol)); 756 ierr = PetscOptionsInt("-pc_hypre_ams_relax_type","Relaxation type for AMS smoother","None",jac->as_relax_type,&jac->as_relax_type,&flag);CHKERRQ(ierr); 757 ierr = PetscOptionsInt("-pc_hypre_ams_relax_times","Number of relaxation steps for AMS smoother","None",jac->as_relax_times,&jac->as_relax_times,&flag2);CHKERRQ(ierr); 758 ierr = PetscOptionsReal("-pc_hypre_ams_relax_weight","Relaxation weight for AMS smoother","None",jac->as_relax_weight,&jac->as_relax_weight,&flag3);CHKERRQ(ierr); 759 ierr = PetscOptionsReal("-pc_hypre_ams_omega","SSOR coefficient for AMS smoother","None",jac->as_omega,&jac->as_omega,&flag4);CHKERRQ(ierr); 760 if (flag || flag2 || flag3 || flag4) { 761 PetscStackCallStandard(HYPRE_AMSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type, 762 jac->as_relax_times, 763 jac->as_relax_weight, 764 jac->as_omega)); 765 } 766 ierr = PetscOptionsReal("-pc_hypre_ams_amg_alpha_theta","Threshold for strong coupling of vector Poisson AMG solver","None",jac->as_amg_alpha_theta,&jac->as_amg_alpha_theta,&flag);CHKERRQ(ierr); 767 n = 5; 768 ierr = PetscOptionsIntArray("-pc_hypre_ams_amg_alpha_options","AMG options for vector Poisson","None",jac->as_amg_alpha_opts,&n,&flag2);CHKERRQ(ierr); 769 if (flag || flag2) { 770 PetscStackCallStandard(HYPRE_AMSSetAlphaAMGOptions,(jac->hsolver,jac->as_amg_alpha_opts[0], /* AMG coarsen type */ 771 jac->as_amg_alpha_opts[1], /* AMG agg_levels */ 772 jac->as_amg_alpha_opts[2], /* AMG relax_type */ 773 jac->as_amg_alpha_theta, 774 jac->as_amg_alpha_opts[3], /* AMG interp_type */ 775 jac->as_amg_alpha_opts[4])); /* AMG Pmax */ 776 } 777 ierr = PetscOptionsReal("-pc_hypre_ams_amg_beta_theta","Threshold for strong coupling of scalar Poisson AMG solver","None",jac->as_amg_beta_theta,&jac->as_amg_beta_theta,&flag);CHKERRQ(ierr); 778 n = 5; 779 ierr = PetscOptionsIntArray("-pc_hypre_ams_amg_beta_options","AMG options for scalar Poisson solver","None",jac->as_amg_beta_opts,&n,&flag2);CHKERRQ(ierr); 780 if (flag || flag2) { 781 PetscStackCallStandard(HYPRE_AMSSetBetaAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */ 782 jac->as_amg_beta_opts[1], /* AMG agg_levels */ 783 jac->as_amg_beta_opts[2], /* AMG relax_type */ 784 jac->as_amg_beta_theta, 785 jac->as_amg_beta_opts[3], /* AMG interp_type */ 786 jac->as_amg_beta_opts[4])); /* AMG Pmax */ 787 } 788 ierr = PetscOptionsTail();CHKERRQ(ierr); 789 PetscFunctionReturn(0); 790 } 791 792 #undef __FUNCT__ 793 #define __FUNCT__ "PCView_HYPRE_AMS" 794 static PetscErrorCode PCView_HYPRE_AMS(PC pc,PetscViewer viewer) 795 { 796 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 797 PetscErrorCode ierr; 798 PetscBool iascii; 799 800 PetscFunctionBegin; 801 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 802 if (iascii) { 803 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS preconditioning\n");CHKERRQ(ierr); 804 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: subspace iterations per application %d\n",jac->as_max_iter);CHKERRQ(ierr); 805 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: subspace cycle type %d\n",jac->ams_cycle_type);CHKERRQ(ierr); 806 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: subspace iteration tolerance %g\n",jac->as_tol);CHKERRQ(ierr); 807 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: smoother type %d\n",jac->as_relax_type);CHKERRQ(ierr); 808 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: number of smoothing steps %d\n",jac->as_relax_times);CHKERRQ(ierr); 809 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: smoother weight %g\n",jac->as_relax_weight);CHKERRQ(ierr); 810 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: smoother omega %g\n",jac->as_omega);CHKERRQ(ierr); 811 if (jac->alpha_Poisson) { 812 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: vector Poisson solver (passed in by user)\n");CHKERRQ(ierr); 813 } else { 814 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: vector Poisson solver (computed) \n");CHKERRQ(ierr); 815 } 816 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG coarsening type %d\n",jac->as_amg_alpha_opts[0]);CHKERRQ(ierr); 817 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG levels of aggressive coarsening %d\n",jac->as_amg_alpha_opts[1]);CHKERRQ(ierr); 818 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG relaxation type %d\n",jac->as_amg_alpha_opts[2]);CHKERRQ(ierr); 819 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG interpolation type %d\n",jac->as_amg_alpha_opts[3]);CHKERRQ(ierr); 820 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG max nonzero elements in interpolation rows %d\n",jac->as_amg_alpha_opts[4]);CHKERRQ(ierr); 821 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG strength threshold %g\n",jac->as_amg_alpha_theta);CHKERRQ(ierr); 822 if (!jac->ams_beta_is_zero) { 823 if (jac->beta_Poisson) { 824 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: scalar Poisson solver (passed in by user)\n");CHKERRQ(ierr); 825 } else { 826 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: scalar Poisson solver (computed) \n");CHKERRQ(ierr); 827 } 828 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG coarsening type %d\n",jac->as_amg_beta_opts[0]);CHKERRQ(ierr); 829 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG levels of aggressive coarsening %d\n",jac->as_amg_beta_opts[1]);CHKERRQ(ierr); 830 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG relaxation type %d\n",jac->as_amg_beta_opts[2]);CHKERRQ(ierr); 831 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG interpolation type %d\n",jac->as_amg_beta_opts[3]);CHKERRQ(ierr); 832 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG max nonzero elements in interpolation rows %d\n",jac->as_amg_beta_opts[4]);CHKERRQ(ierr); 833 ierr = PetscViewerASCIIPrintf(viewer," HYPRE AMS: boomerAMG strength threshold %g\n",jac->as_amg_beta_theta);CHKERRQ(ierr); 834 } 835 } 836 PetscFunctionReturn(0); 837 } 838 839 #undef __FUNCT__ 840 #define __FUNCT__ "PCSetFromOptions_HYPRE_ADS" 841 static PetscErrorCode PCSetFromOptions_HYPRE_ADS(PetscOptions *PetscOptionsObject,PC pc) 842 { 843 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 844 PetscErrorCode ierr; 845 PetscInt n; 846 PetscBool flag,flag2,flag3,flag4; 847 848 PetscFunctionBegin; 849 ierr = PetscOptionsHead(PetscOptionsObject,"HYPRE ADS Options");CHKERRQ(ierr); 850 ierr = PetscOptionsInt("-pc_hypre_ads_print_level","Debugging output level for ADS","None",jac->as_print,&jac->as_print,&flag);CHKERRQ(ierr); 851 if (flag) PetscStackCallStandard(HYPRE_ADSSetPrintLevel,(jac->hsolver,jac->as_print)); 852 ierr = PetscOptionsInt("-pc_hypre_ads_max_iter","Maximum number of ADS multigrid iterations within PCApply","None",jac->as_max_iter,&jac->as_max_iter,&flag);CHKERRQ(ierr); 853 if (flag) PetscStackCallStandard(HYPRE_ADSSetMaxIter,(jac->hsolver,jac->as_max_iter)); 854 ierr = PetscOptionsInt("-pc_hypre_ads_cycle_type","Cycle type for ADS multigrid","None",jac->ads_cycle_type,&jac->ads_cycle_type,&flag);CHKERRQ(ierr); 855 if (flag) PetscStackCallStandard(HYPRE_ADSSetCycleType,(jac->hsolver,jac->ads_cycle_type)); 856 ierr = PetscOptionsReal("-pc_hypre_ads_tol","Error tolerance for ADS multigrid","None",jac->as_tol,&jac->as_tol,&flag);CHKERRQ(ierr); 857 if (flag) PetscStackCallStandard(HYPRE_ADSSetTol,(jac->hsolver,jac->as_tol)); 858 ierr = PetscOptionsInt("-pc_hypre_ads_relax_type","Relaxation type for ADS smoother","None",jac->as_relax_type,&jac->as_relax_type,&flag);CHKERRQ(ierr); 859 ierr = PetscOptionsInt("-pc_hypre_ads_relax_times","Number of relaxation steps for ADS smoother","None",jac->as_relax_times,&jac->as_relax_times,&flag2);CHKERRQ(ierr); 860 ierr = PetscOptionsReal("-pc_hypre_ads_relax_weight","Relaxation weight for ADS smoother","None",jac->as_relax_weight,&jac->as_relax_weight,&flag3);CHKERRQ(ierr); 861 ierr = PetscOptionsReal("-pc_hypre_ads_omega","SSOR coefficient for ADS smoother","None",jac->as_omega,&jac->as_omega,&flag4);CHKERRQ(ierr); 862 if (flag || flag2 || flag3 || flag4) { 863 PetscStackCallStandard(HYPRE_ADSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type, 864 jac->as_relax_times, 865 jac->as_relax_weight, 866 jac->as_omega)); 867 } 868 ierr = PetscOptionsReal("-pc_hypre_ads_ams_theta","Threshold for strong coupling of AMS solver inside ADS","None",jac->as_amg_alpha_theta,&jac->as_amg_alpha_theta,&flag);CHKERRQ(ierr); 869 n = 5; 870 ierr = PetscOptionsIntArray("-pc_hypre_ads_ams_options","AMG options for AMS solver inside ADS","None",jac->as_amg_alpha_opts,&n,&flag2);CHKERRQ(ierr); 871 ierr = PetscOptionsInt("-pc_hypre_ads_ams_cycle_type","Cycle type for AMS solver inside ADS","None",jac->ams_cycle_type,&jac->ams_cycle_type,&flag3);CHKERRQ(ierr); 872 if (flag || flag2 || flag3) { 873 PetscStackCallStandard(HYPRE_ADSSetAMSOptions,(jac->hsolver,jac->ams_cycle_type, /* AMS cycle type */ 874 jac->as_amg_alpha_opts[0], /* AMG coarsen type */ 875 jac->as_amg_alpha_opts[1], /* AMG agg_levels */ 876 jac->as_amg_alpha_opts[2], /* AMG relax_type */ 877 jac->as_amg_alpha_theta, 878 jac->as_amg_alpha_opts[3], /* AMG interp_type */ 879 jac->as_amg_alpha_opts[4])); /* AMG Pmax */ 880 } 881 ierr = PetscOptionsReal("-pc_hypre_ads_amg_theta","Threshold for strong coupling of vector AMG solver inside ADS","None",jac->as_amg_beta_theta,&jac->as_amg_beta_theta,&flag);CHKERRQ(ierr); 882 n = 5; 883 ierr = PetscOptionsIntArray("-pc_hypre_ads_amg_options","AMG options for vector AMG solver inside ADS","None",jac->as_amg_beta_opts,&n,&flag2);CHKERRQ(ierr); 884 if (flag || flag2) { 885 PetscStackCallStandard(HYPRE_ADSSetAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */ 886 jac->as_amg_beta_opts[1], /* AMG agg_levels */ 887 jac->as_amg_beta_opts[2], /* AMG relax_type */ 888 jac->as_amg_beta_theta, 889 jac->as_amg_beta_opts[3], /* AMG interp_type */ 890 jac->as_amg_beta_opts[4])); /* AMG Pmax */ 891 } 892 ierr = PetscOptionsTail();CHKERRQ(ierr); 893 PetscFunctionReturn(0); 894 } 895 896 #undef __FUNCT__ 897 #define __FUNCT__ "PCView_HYPRE_ADS" 898 static PetscErrorCode PCView_HYPRE_ADS(PC pc,PetscViewer viewer) 899 { 900 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 901 PetscErrorCode ierr; 902 PetscBool iascii; 903 904 PetscFunctionBegin; 905 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 906 if (iascii) { 907 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS preconditioning\n");CHKERRQ(ierr); 908 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: subspace iterations per application %d\n",jac->as_max_iter);CHKERRQ(ierr); 909 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: subspace cycle type %d\n",jac->ads_cycle_type);CHKERRQ(ierr); 910 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: subspace iteration tolerance %g\n",jac->as_tol);CHKERRQ(ierr); 911 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: smoother type %d\n",jac->as_relax_type);CHKERRQ(ierr); 912 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: number of smoothing steps %d\n",jac->as_relax_times);CHKERRQ(ierr); 913 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: smoother weight %g\n",jac->as_relax_weight);CHKERRQ(ierr); 914 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: smoother omega %g\n",jac->as_omega);CHKERRQ(ierr); 915 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: AMS solver\n");CHKERRQ(ierr); 916 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: subspace cycle type %d\n",jac->ams_cycle_type);CHKERRQ(ierr); 917 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG coarsening type %d\n",jac->as_amg_alpha_opts[0]);CHKERRQ(ierr); 918 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG levels of aggressive coarsening %d\n",jac->as_amg_alpha_opts[1]);CHKERRQ(ierr); 919 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG relaxation type %d\n",jac->as_amg_alpha_opts[2]);CHKERRQ(ierr); 920 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG interpolation type %d\n",jac->as_amg_alpha_opts[3]);CHKERRQ(ierr); 921 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG max nonzero elements in interpolation rows %d\n",jac->as_amg_alpha_opts[4]);CHKERRQ(ierr); 922 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG strength threshold %g\n",jac->as_amg_alpha_theta);CHKERRQ(ierr); 923 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: vector Poisson solver\n");CHKERRQ(ierr); 924 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG coarsening type %d\n",jac->as_amg_beta_opts[0]);CHKERRQ(ierr); 925 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG levels of aggressive coarsening %d\n",jac->as_amg_beta_opts[1]);CHKERRQ(ierr); 926 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG relaxation type %d\n",jac->as_amg_beta_opts[2]);CHKERRQ(ierr); 927 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG interpolation type %d\n",jac->as_amg_beta_opts[3]);CHKERRQ(ierr); 928 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG max nonzero elements in interpolation rows %d\n",jac->as_amg_beta_opts[4]);CHKERRQ(ierr); 929 ierr = PetscViewerASCIIPrintf(viewer," HYPRE ADS: boomerAMG strength threshold %g\n",jac->as_amg_beta_theta);CHKERRQ(ierr); 930 } 931 PetscFunctionReturn(0); 932 } 933 934 #undef __FUNCT__ 935 #define __FUNCT__ "PCHYPRESetDiscreteGradient_HYPRE" 936 static PetscErrorCode PCHYPRESetDiscreteGradient_HYPRE(PC pc, Mat G) 937 { 938 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 939 HYPRE_ParCSRMatrix parcsr_G; 940 PetscErrorCode ierr; 941 942 PetscFunctionBegin; 943 /* throw away any discrete gradient if already set */ 944 if (jac->G) PetscStackCallStandard(HYPRE_IJMatrixDestroy,(jac->G)); 945 ierr = MatHYPRE_IJMatrixCreate(G,&jac->G);CHKERRQ(ierr); 946 ierr = MatHYPRE_IJMatrixCopy(G,jac->G);CHKERRQ(ierr); 947 PetscStackCallStandard(HYPRE_IJMatrixGetObject,(jac->G,(void**)(&parcsr_G))); 948 PetscStackCall("Hypre set gradient",ierr = (*jac->setdgrad)(jac->hsolver,parcsr_G);CHKERRQ(ierr);); 949 PetscFunctionReturn(0); 950 } 951 952 #undef __FUNCT__ 953 #define __FUNCT__ "PCHYPRESetDiscreteGradient" 954 /*@ 955 PCHYPRESetDiscreteGradient - Set discrete gradient matrix 956 957 Collective on PC 958 959 Input Parameters: 960 + pc - the preconditioning context 961 - G - the discrete gradient 962 963 Level: intermediate 964 965 Notes: G should have as many rows as the number of edges and as many columns as the number of vertices in the mesh 966 Each row of G has 2 nonzeros, with column indexes being the global indexes of edge's endpoints: matrix entries are +1 and -1 depending on edge orientation 967 968 .seealso: 969 @*/ 970 PetscErrorCode PCHYPRESetDiscreteGradient(PC pc, Mat G) 971 { 972 PetscErrorCode ierr; 973 974 PetscFunctionBegin; 975 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 976 PetscValidHeaderSpecific(G,MAT_CLASSID,2); 977 PetscCheckSameComm(pc,1,G,2); 978 ierr = PetscTryMethod(pc,"PCHYPRESetDiscreteGradient_C",(PC,Mat),(pc,G));CHKERRQ(ierr); 979 PetscFunctionReturn(0); 980 } 981 982 #undef __FUNCT__ 983 #define __FUNCT__ "PCHYPRESetDiscreteCurl_HYPRE" 984 static PetscErrorCode PCHYPRESetDiscreteCurl_HYPRE(PC pc, Mat C) 985 { 986 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 987 HYPRE_ParCSRMatrix parcsr_C; 988 PetscErrorCode ierr; 989 990 PetscFunctionBegin; 991 /* throw away any discrete curl if already set */ 992 if (jac->C) PetscStackCallStandard(HYPRE_IJMatrixDestroy,(jac->C)); 993 ierr = MatHYPRE_IJMatrixCreate(C,&jac->C);CHKERRQ(ierr); 994 ierr = MatHYPRE_IJMatrixCopy(C,jac->C);CHKERRQ(ierr); 995 PetscStackCallStandard(HYPRE_IJMatrixGetObject,(jac->C,(void**)(&parcsr_C))); 996 PetscStackCall("Hypre set curl",ierr = (*jac->setdcurl)(jac->hsolver,parcsr_C);CHKERRQ(ierr);); 997 PetscFunctionReturn(0); 998 } 999 1000 #undef __FUNCT__ 1001 #define __FUNCT__ "PCHYPRESetDiscreteCurl" 1002 /*@ 1003 PCHYPRESetDiscreteCurl - Set discrete curl matrix 1004 1005 Collective on PC 1006 1007 Input Parameters: 1008 + pc - the preconditioning context 1009 - C - the discrete curl 1010 1011 Level: intermediate 1012 1013 Notes: C should have as many rows as the number of faces and as many columns as the number of edges in the mesh 1014 Each row of G has as many nonzeros as the number of edges of a face, with column indexes being the global indexes of the corresponding edge: matrix entries are +1 and -1 depending on edge orientation with respect to the face orientation 1015 1016 .seealso: 1017 @*/ 1018 PetscErrorCode PCHYPRESetDiscreteCurl(PC pc, Mat C) 1019 { 1020 PetscErrorCode ierr; 1021 1022 PetscFunctionBegin; 1023 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1024 PetscValidHeaderSpecific(C,MAT_CLASSID,2); 1025 PetscCheckSameComm(pc,1,C,2); 1026 ierr = PetscTryMethod(pc,"PCHYPRESetDiscreteCurl_C",(PC,Mat),(pc,C));CHKERRQ(ierr); 1027 PetscFunctionReturn(0); 1028 } 1029 1030 #undef __FUNCT__ 1031 #define __FUNCT__ "PCHYPRESetAlphaPoissonMatrix_HYPRE_AMS" 1032 static PetscErrorCode PCHYPRESetAlphaPoissonMatrix_HYPRE_AMS(PC pc, Mat A) 1033 { 1034 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1035 HYPRE_ParCSRMatrix parcsr_alpha_Poisson; 1036 PetscErrorCode ierr; 1037 1038 PetscFunctionBegin; 1039 /* throw away any matrix if already set */ 1040 if (jac->alpha_Poisson) PetscStackCallStandard(HYPRE_IJMatrixDestroy,(jac->alpha_Poisson)); 1041 ierr = MatHYPRE_IJMatrixCreate(A,&jac->alpha_Poisson);CHKERRQ(ierr); 1042 ierr = MatHYPRE_IJMatrixCopy(A,jac->alpha_Poisson);CHKERRQ(ierr); 1043 PetscStackCallStandard(HYPRE_IJMatrixGetObject,(jac->alpha_Poisson,(void**)(&parcsr_alpha_Poisson))); 1044 PetscStackCallStandard(HYPRE_AMSSetAlphaPoissonMatrix,(jac->hsolver,parcsr_alpha_Poisson)); 1045 PetscFunctionReturn(0); 1046 } 1047 1048 #undef __FUNCT__ 1049 #define __FUNCT__ "PCHYPRESetAlphaPoissonMatrix" 1050 /*@ 1051 PCHYPRESetAlphaPoissonMatrix - Set vector Poisson matrix 1052 1053 Collective on PC 1054 1055 Input Parameters: 1056 + pc - the preconditioning context 1057 - A - the matrix 1058 1059 Level: intermediate 1060 1061 Notes: A should be obtained by discretizing the vector valued Poisson problem with linear finite elements 1062 1063 .seealso: 1064 @*/ 1065 PetscErrorCode PCHYPRESetAlphaPoissonMatrix(PC pc, Mat A) 1066 { 1067 PetscErrorCode ierr; 1068 1069 PetscFunctionBegin; 1070 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1071 PetscValidHeaderSpecific(A,MAT_CLASSID,2); 1072 PetscCheckSameComm(pc,1,A,2); 1073 ierr = PetscTryMethod(pc,"PCHYPRESetAlphaPoissonMatrix_C",(PC,Mat),(pc,A));CHKERRQ(ierr); 1074 PetscFunctionReturn(0); 1075 } 1076 1077 #undef __FUNCT__ 1078 #define __FUNCT__ "PCHYPRESetBetaPoissonMatrix_HYPRE_AMS" 1079 static PetscErrorCode PCHYPRESetBetaPoissonMatrix_HYPRE_AMS(PC pc, Mat A) 1080 { 1081 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1082 HYPRE_ParCSRMatrix parcsr_beta_Poisson; 1083 PetscErrorCode ierr; 1084 1085 PetscFunctionBegin; 1086 if (!A) { 1087 PetscStackCallStandard(HYPRE_AMSSetBetaPoissonMatrix,(jac->hsolver,NULL)); 1088 jac->ams_beta_is_zero = PETSC_TRUE; 1089 PetscFunctionReturn(0); 1090 } 1091 jac->ams_beta_is_zero = PETSC_FALSE; 1092 /* throw away any matrix if already set */ 1093 if (jac->beta_Poisson) PetscStackCallStandard(HYPRE_IJMatrixDestroy,(jac->beta_Poisson)); 1094 ierr = MatHYPRE_IJMatrixCreate(A,&jac->beta_Poisson);CHKERRQ(ierr); 1095 ierr = MatHYPRE_IJMatrixCopy(A,jac->beta_Poisson);CHKERRQ(ierr); 1096 PetscStackCallStandard(HYPRE_IJMatrixGetObject,(jac->beta_Poisson,(void**)(&parcsr_beta_Poisson))); 1097 PetscStackCallStandard(HYPRE_AMSSetBetaPoissonMatrix,(jac->hsolver,parcsr_beta_Poisson)); 1098 PetscFunctionReturn(0); 1099 } 1100 1101 #undef __FUNCT__ 1102 #define __FUNCT__ "PCHYPRESetBetaPoissonMatrix" 1103 /*@ 1104 PCHYPRESetBetaPoissonMatrix - Set Poisson matrix 1105 1106 Collective on PC 1107 1108 Input Parameters: 1109 + pc - the preconditioning context 1110 - A - the matrix 1111 1112 Level: intermediate 1113 1114 Notes: A should be obtained by discretizing the Poisson problem with linear finite elements. 1115 Following HYPRE convention, the scalar Poisson solver of AMS can be turned off by passing NULL. 1116 1117 .seealso: 1118 @*/ 1119 PetscErrorCode PCHYPRESetBetaPoissonMatrix(PC pc, Mat A) 1120 { 1121 PetscErrorCode ierr; 1122 1123 PetscFunctionBegin; 1124 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1125 if (A) { 1126 PetscValidHeaderSpecific(A,MAT_CLASSID,2); 1127 PetscCheckSameComm(pc,1,A,2); 1128 } 1129 ierr = PetscTryMethod(pc,"PCHYPRESetBetaPoissonMatrix_C",(PC,Mat),(pc,A));CHKERRQ(ierr); 1130 PetscFunctionReturn(0); 1131 } 1132 1133 #undef __FUNCT__ 1134 #define __FUNCT__ "PCHYPRESetEdgeConstantVectors_HYPRE_AMS" 1135 static PetscErrorCode PCHYPRESetEdgeConstantVectors_HYPRE_AMS(PC pc,Vec ozz, Vec zoz, Vec zzo) 1136 { 1137 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1138 HYPRE_ParVector par_ozz,par_zoz,par_zzo; 1139 PetscInt dim; 1140 PetscErrorCode ierr; 1141 1142 PetscFunctionBegin; 1143 /* throw away any vector if already set */ 1144 if (jac->constants[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[0])); 1145 if (jac->constants[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[1])); 1146 if (jac->constants[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->constants[2])); 1147 jac->constants[0] = NULL; 1148 jac->constants[1] = NULL; 1149 jac->constants[2] = NULL; 1150 ierr = VecHYPRE_IJVectorCreate(ozz,&jac->constants[0]);CHKERRQ(ierr); 1151 ierr = VecHYPRE_IJVectorCopy(ozz,jac->constants[0]);CHKERRQ(ierr); 1152 PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->constants[0],(void**)(&par_ozz))); 1153 ierr = VecHYPRE_IJVectorCreate(zoz,&jac->constants[1]);CHKERRQ(ierr); 1154 ierr = VecHYPRE_IJVectorCopy(zoz,jac->constants[1]);CHKERRQ(ierr); 1155 PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->constants[1],(void**)(&par_zoz))); 1156 dim = 2; 1157 if (zzo) { 1158 ierr = VecHYPRE_IJVectorCreate(zzo,&jac->constants[2]);CHKERRQ(ierr); 1159 ierr = VecHYPRE_IJVectorCopy(zzo,jac->constants[2]);CHKERRQ(ierr); 1160 PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->constants[2],(void**)(&par_zzo))); 1161 dim++; 1162 } 1163 PetscStackCallStandard(HYPRE_AMSSetEdgeConstantVectors,(jac->hsolver,par_ozz,par_zoz,par_zzo)); 1164 PetscStackCallStandard(HYPRE_AMSSetDimension,(jac->hsolver,dim)); 1165 PetscFunctionReturn(0); 1166 } 1167 1168 #undef __FUNCT__ 1169 #define __FUNCT__ "PCHYPRESetEdgeConstantVectors" 1170 /*@ 1171 PCHYPRESetEdgeConstantVectors - Set the representation of the constant vector fields in edge element basis 1172 1173 Collective on PC 1174 1175 Input Parameters: 1176 + pc - the preconditioning context 1177 - ozz - vector representing (1,0,0) (or (1,0) in 2D) 1178 - zoz - vector representing (0,1,0) (or (0,1) in 2D) 1179 - zzo - vector representing (0,0,1) (use NULL in 2D) 1180 1181 Level: intermediate 1182 1183 Notes: 1184 1185 .seealso: 1186 @*/ 1187 PetscErrorCode PCHYPRESetEdgeConstantVectors(PC pc, Vec ozz, Vec zoz, Vec zzo) 1188 { 1189 PetscErrorCode ierr; 1190 1191 PetscFunctionBegin; 1192 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1193 PetscValidHeaderSpecific(ozz,VEC_CLASSID,2); 1194 PetscValidHeaderSpecific(zoz,VEC_CLASSID,3); 1195 if (zzo) PetscValidHeaderSpecific(zzo,VEC_CLASSID,4); 1196 PetscCheckSameComm(pc,1,ozz,2); 1197 PetscCheckSameComm(pc,1,zoz,3); 1198 if (zzo) PetscCheckSameComm(pc,1,zzo,4); 1199 ierr = PetscTryMethod(pc,"PCHYPRESetEdgeConstantVectors_C",(PC,Vec,Vec,Vec),(pc,ozz,zoz,zzo));CHKERRQ(ierr); 1200 PetscFunctionReturn(0); 1201 } 1202 1203 #undef __FUNCT__ 1204 #define __FUNCT__ "PCSetCoordinates_HYPRE" 1205 static PetscErrorCode PCSetCoordinates_HYPRE(PC pc, PetscInt dim, PetscInt nloc, PetscReal *coords) 1206 { 1207 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1208 Vec tv; 1209 HYPRE_ParVector par_coords[3]; 1210 PetscInt i; 1211 PetscErrorCode ierr; 1212 1213 PetscFunctionBegin; 1214 /* throw away any coordinate vector if already set */ 1215 if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[0])); 1216 if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[1])); 1217 if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorDestroy,(jac->coords[2])); 1218 /* set problem's dimension */ 1219 if (jac->setdim) { 1220 PetscStackCall("Hypre set dim",ierr = (*jac->setdim)(jac->hsolver,dim);CHKERRQ(ierr);); 1221 } 1222 /* compute IJ vector for coordinates */ 1223 ierr = VecCreate(PetscObjectComm((PetscObject)pc),&tv);CHKERRQ(ierr); 1224 ierr = VecSetType(tv,VECSTANDARD);CHKERRQ(ierr); 1225 ierr = VecSetSizes(tv,nloc,PETSC_DECIDE);CHKERRQ(ierr); 1226 for (i=0;i<dim;i++) { 1227 PetscScalar *array; 1228 PetscInt j; 1229 1230 ierr = VecHYPRE_IJVectorCreate(tv,&jac->coords[i]);CHKERRQ(ierr); 1231 ierr = VecGetArray(tv,&array);CHKERRQ(ierr); 1232 for (j=0;j<nloc;j++) { 1233 array[j] = coords[j*dim+i]; 1234 } 1235 PetscStackCallStandard(HYPRE_IJVectorSetValues,(jac->coords[i],nloc,NULL,array)); 1236 PetscStackCallStandard(HYPRE_IJVectorAssemble,(jac->coords[i])); 1237 ierr = VecRestoreArray(tv,&array);CHKERRQ(ierr); 1238 } 1239 ierr = VecDestroy(&tv);CHKERRQ(ierr); 1240 /* pass parCSR vectors to AMS solver */ 1241 par_coords[0] = NULL; 1242 par_coords[1] = NULL; 1243 par_coords[2] = NULL; 1244 if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[0],(void**)(&par_coords[0]))); 1245 if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[1],(void**)(&par_coords[1]))); 1246 if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorGetObject,(jac->coords[2],(void**)(&par_coords[2]))); 1247 PetscStackCall("Hypre set coords",ierr = (*jac->setcoord)(jac->hsolver,par_coords[0],par_coords[1],par_coords[2]);CHKERRQ(ierr);); 1248 PetscFunctionReturn(0); 1249 } 1250 1251 /* ---------------------------------------------------------------------------------*/ 1252 1253 #undef __FUNCT__ 1254 #define __FUNCT__ "PCHYPREGetType_HYPRE" 1255 static PetscErrorCode PCHYPREGetType_HYPRE(PC pc,const char *name[]) 1256 { 1257 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1258 1259 PetscFunctionBegin; 1260 *name = jac->hypre_type; 1261 PetscFunctionReturn(0); 1262 } 1263 1264 #undef __FUNCT__ 1265 #define __FUNCT__ "PCHYPRESetType_HYPRE" 1266 static PetscErrorCode PCHYPRESetType_HYPRE(PC pc,const char name[]) 1267 { 1268 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1269 PetscErrorCode ierr; 1270 PetscBool flag; 1271 1272 PetscFunctionBegin; 1273 if (jac->hypre_type) { 1274 ierr = PetscStrcmp(jac->hypre_type,name,&flag);CHKERRQ(ierr); 1275 if (!flag) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ORDER,"Cannot reset the HYPRE preconditioner type once it has been set"); 1276 PetscFunctionReturn(0); 1277 } else { 1278 ierr = PetscStrallocpy(name, &jac->hypre_type);CHKERRQ(ierr); 1279 } 1280 1281 jac->maxiter = PETSC_DEFAULT; 1282 jac->tol = PETSC_DEFAULT; 1283 jac->printstatistics = PetscLogPrintInfo; 1284 1285 ierr = PetscStrcmp("pilut",jac->hypre_type,&flag);CHKERRQ(ierr); 1286 if (flag) { 1287 PetscStackCallStandard(HYPRE_ParCSRPilutCreate,(jac->comm_hypre,&jac->hsolver)); 1288 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Pilut; 1289 pc->ops->view = PCView_HYPRE_Pilut; 1290 jac->destroy = HYPRE_ParCSRPilutDestroy; 1291 jac->setup = HYPRE_ParCSRPilutSetup; 1292 jac->solve = HYPRE_ParCSRPilutSolve; 1293 jac->factorrowsize = PETSC_DEFAULT; 1294 PetscFunctionReturn(0); 1295 } 1296 ierr = PetscStrcmp("parasails",jac->hypre_type,&flag);CHKERRQ(ierr); 1297 if (flag) { 1298 PetscStackCallStandard(HYPRE_ParaSailsCreate,(jac->comm_hypre,&jac->hsolver)); 1299 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ParaSails; 1300 pc->ops->view = PCView_HYPRE_ParaSails; 1301 jac->destroy = HYPRE_ParaSailsDestroy; 1302 jac->setup = HYPRE_ParaSailsSetup; 1303 jac->solve = HYPRE_ParaSailsSolve; 1304 /* initialize */ 1305 jac->nlevels = 1; 1306 jac->threshhold = .1; 1307 jac->filter = .1; 1308 jac->loadbal = 0; 1309 if (PetscLogPrintInfo) jac->logging = (int) PETSC_TRUE; 1310 else jac->logging = (int) PETSC_FALSE; 1311 1312 jac->ruse = (int) PETSC_FALSE; 1313 jac->symt = 0; 1314 PetscStackCallStandard(HYPRE_ParaSailsSetParams,(jac->hsolver,jac->threshhold,jac->nlevels)); 1315 PetscStackCallStandard(HYPRE_ParaSailsSetFilter,(jac->hsolver,jac->filter)); 1316 PetscStackCallStandard(HYPRE_ParaSailsSetLoadbal,(jac->hsolver,jac->loadbal)); 1317 PetscStackCallStandard(HYPRE_ParaSailsSetLogging,(jac->hsolver,jac->logging)); 1318 PetscStackCallStandard(HYPRE_ParaSailsSetReuse,(jac->hsolver,jac->ruse)); 1319 PetscStackCallStandard(HYPRE_ParaSailsSetSym,(jac->hsolver,jac->symt)); 1320 PetscFunctionReturn(0); 1321 } 1322 ierr = PetscStrcmp("boomeramg",jac->hypre_type,&flag);CHKERRQ(ierr); 1323 if (flag) { 1324 ierr = HYPRE_BoomerAMGCreate(&jac->hsolver); 1325 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_BoomerAMG; 1326 pc->ops->view = PCView_HYPRE_BoomerAMG; 1327 pc->ops->applytranspose = PCApplyTranspose_HYPRE_BoomerAMG; 1328 pc->ops->applyrichardson = PCApplyRichardson_HYPRE_BoomerAMG; 1329 jac->destroy = HYPRE_BoomerAMGDestroy; 1330 jac->setup = HYPRE_BoomerAMGSetup; 1331 jac->solve = HYPRE_BoomerAMGSolve; 1332 jac->applyrichardson = PETSC_FALSE; 1333 /* these defaults match the hypre defaults */ 1334 jac->cycletype = 1; 1335 jac->maxlevels = 25; 1336 jac->maxiter = 1; 1337 jac->tol = 0.0; /* tolerance of zero indicates use as preconditioner (suppresses convergence errors) */ 1338 jac->truncfactor = 0.0; 1339 jac->strongthreshold = .25; 1340 jac->maxrowsum = .9; 1341 jac->coarsentype = 6; 1342 jac->measuretype = 0; 1343 jac->gridsweeps[0] = jac->gridsweeps[1] = jac->gridsweeps[2] = 1; 1344 jac->relaxtype[0] = jac->relaxtype[1] = 6; /* Defaults to SYMMETRIC since in PETSc we are using a a PC - most likely with CG */ 1345 jac->relaxtype[2] = 9; /*G.E. */ 1346 jac->relaxweight = 1.0; 1347 jac->outerrelaxweight = 1.0; 1348 jac->relaxorder = 1; 1349 jac->interptype = 0; 1350 jac->agg_nl = 0; 1351 jac->pmax = 0; 1352 jac->truncfactor = 0.0; 1353 jac->agg_num_paths = 1; 1354 1355 jac->nodal_coarsen = 0; 1356 jac->nodal_relax = PETSC_FALSE; 1357 jac->nodal_relax_levels = 1; 1358 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleType,(jac->hsolver,jac->cycletype)); 1359 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxLevels,(jac->hsolver,jac->maxlevels)); 1360 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,(jac->hsolver,jac->maxiter)); 1361 PetscStackCallStandard(HYPRE_BoomerAMGSetTol,(jac->hsolver,jac->tol)); 1362 PetscStackCallStandard(HYPRE_BoomerAMGSetTruncFactor,(jac->hsolver,jac->truncfactor)); 1363 PetscStackCallStandard(HYPRE_BoomerAMGSetStrongThreshold,(jac->hsolver,jac->strongthreshold)); 1364 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxRowSum,(jac->hsolver,jac->maxrowsum)); 1365 PetscStackCallStandard(HYPRE_BoomerAMGSetCoarsenType,(jac->hsolver,jac->coarsentype)); 1366 PetscStackCallStandard(HYPRE_BoomerAMGSetMeasureType,(jac->hsolver,jac->measuretype)); 1367 PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxOrder,(jac->hsolver, jac->relaxorder)); 1368 PetscStackCallStandard(HYPRE_BoomerAMGSetInterpType,(jac->hsolver,jac->interptype)); 1369 PetscStackCallStandard(HYPRE_BoomerAMGSetAggNumLevels,(jac->hsolver,jac->agg_nl)); 1370 PetscStackCallStandard(HYPRE_BoomerAMGSetPMaxElmts,(jac->hsolver,jac->pmax)); 1371 PetscStackCallStandard(HYPRE_BoomerAMGSetNumPaths,(jac->hsolver,jac->agg_num_paths)); 1372 PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxType,(jac->hsolver, jac->relaxtype[0])); /*defaults coarse to 9*/ 1373 PetscStackCallStandard(HYPRE_BoomerAMGSetNumSweeps,(jac->hsolver, jac->gridsweeps[0])); /*defaults coarse to 1 */ 1374 PetscFunctionReturn(0); 1375 } 1376 ierr = PetscStrcmp("ams",jac->hypre_type,&flag);CHKERRQ(ierr); 1377 if (flag) { 1378 ierr = HYPRE_AMSCreate(&jac->hsolver); 1379 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_AMS; 1380 pc->ops->view = PCView_HYPRE_AMS; 1381 jac->destroy = HYPRE_AMSDestroy; 1382 jac->setup = HYPRE_AMSSetup; 1383 jac->solve = HYPRE_AMSSolve; 1384 jac->setdgrad = HYPRE_AMSSetDiscreteGradient; 1385 jac->setcoord = HYPRE_AMSSetCoordinateVectors; 1386 jac->setdim = HYPRE_AMSSetDimension; 1387 jac->coords[0] = NULL; 1388 jac->coords[1] = NULL; 1389 jac->coords[2] = NULL; 1390 jac->G = NULL; 1391 /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE AMS */ 1392 jac->as_print = 0; 1393 jac->as_max_iter = 1; /* used as a preconditioner */ 1394 jac->as_tol = 0.; /* used as a preconditioner */ 1395 jac->ams_cycle_type = 13; 1396 /* Smoothing options */ 1397 jac->as_relax_type = 2; 1398 jac->as_relax_times = 1; 1399 jac->as_relax_weight = 1.0; 1400 jac->as_omega = 1.0; 1401 /* Vector valued Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */ 1402 jac->as_amg_alpha_opts[0] = 10; 1403 jac->as_amg_alpha_opts[1] = 1; 1404 jac->as_amg_alpha_opts[2] = 8; 1405 jac->as_amg_alpha_opts[3] = 6; 1406 jac->as_amg_alpha_opts[4] = 4; 1407 jac->as_amg_alpha_theta = 0.25; 1408 /* Scalar Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */ 1409 jac->ams_beta_is_zero = PETSC_FALSE; 1410 jac->as_amg_beta_opts[0] = 10; 1411 jac->as_amg_beta_opts[1] = 1; 1412 jac->as_amg_beta_opts[2] = 8; 1413 jac->as_amg_beta_opts[3] = 6; 1414 jac->as_amg_beta_opts[4] = 4; 1415 jac->as_amg_beta_theta = 0.25; 1416 PetscStackCallStandard(HYPRE_AMSSetPrintLevel,(jac->hsolver,jac->as_print)); 1417 PetscStackCallStandard(HYPRE_AMSSetMaxIter,(jac->hsolver,jac->as_max_iter)); 1418 PetscStackCallStandard(HYPRE_AMSSetCycleType,(jac->hsolver,jac->ams_cycle_type)); 1419 PetscStackCallStandard(HYPRE_AMSSetTol,(jac->hsolver,jac->as_tol)); 1420 PetscStackCallStandard(HYPRE_AMSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type, 1421 jac->as_relax_times, 1422 jac->as_relax_weight, 1423 jac->as_omega)); 1424 PetscStackCallStandard(HYPRE_AMSSetAlphaAMGOptions,(jac->hsolver,jac->as_amg_alpha_opts[0], /* AMG coarsen type */ 1425 jac->as_amg_alpha_opts[1], /* AMG agg_levels */ 1426 jac->as_amg_alpha_opts[2], /* AMG relax_type */ 1427 jac->as_amg_alpha_theta, 1428 jac->as_amg_alpha_opts[3], /* AMG interp_type */ 1429 jac->as_amg_alpha_opts[4])); /* AMG Pmax */ 1430 PetscStackCallStandard(HYPRE_AMSSetBetaAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */ 1431 jac->as_amg_beta_opts[1], /* AMG agg_levels */ 1432 jac->as_amg_beta_opts[2], /* AMG relax_type */ 1433 jac->as_amg_beta_theta, 1434 jac->as_amg_beta_opts[3], /* AMG interp_type */ 1435 jac->as_amg_beta_opts[4])); /* AMG Pmax */ 1436 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCSetCoordinates_C",PCSetCoordinates_HYPRE);CHKERRQ(ierr); 1437 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteGradient_C",PCHYPRESetDiscreteGradient_HYPRE);CHKERRQ(ierr); 1438 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetEdgeConstantVectors_C",PCHYPRESetEdgeConstantVectors_HYPRE_AMS);CHKERRQ(ierr); 1439 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetAlphaPoissonMatrix_C",PCHYPRESetAlphaPoissonMatrix_HYPRE_AMS);CHKERRQ(ierr); 1440 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetBetaPoissonMatrix_C",PCHYPRESetBetaPoissonMatrix_HYPRE_AMS);CHKERRQ(ierr); 1441 PetscFunctionReturn(0); 1442 } 1443 ierr = PetscStrcmp("ads",jac->hypre_type,&flag);CHKERRQ(ierr); 1444 if (flag) { 1445 ierr = HYPRE_ADSCreate(&jac->hsolver); 1446 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ADS; 1447 pc->ops->view = PCView_HYPRE_ADS; 1448 jac->destroy = HYPRE_ADSDestroy; 1449 jac->setup = HYPRE_ADSSetup; 1450 jac->solve = HYPRE_ADSSolve; 1451 jac->setdgrad = HYPRE_ADSSetDiscreteGradient; 1452 jac->setdcurl = HYPRE_ADSSetDiscreteCurl; 1453 jac->setcoord = HYPRE_ADSSetCoordinateVectors; 1454 jac->coords[0] = NULL; 1455 jac->coords[1] = NULL; 1456 jac->coords[2] = NULL; 1457 jac->G = NULL; 1458 jac->C = NULL; 1459 /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE ADS */ 1460 jac->as_print = 0; 1461 jac->as_max_iter = 1; /* used as a preconditioner */ 1462 jac->as_tol = 0.; /* used as a preconditioner */ 1463 jac->ads_cycle_type = 13; 1464 /* Smoothing options */ 1465 jac->as_relax_type = 2; 1466 jac->as_relax_times = 1; 1467 jac->as_relax_weight = 1.0; 1468 jac->as_omega = 1.0; 1469 /* AMS solver parameters: cycle_type, coarsen type, agg_levels, relax_type, interp_type, Pmax */ 1470 jac->ams_cycle_type = 14; 1471 jac->as_amg_alpha_opts[0] = 10; 1472 jac->as_amg_alpha_opts[1] = 1; 1473 jac->as_amg_alpha_opts[2] = 6; 1474 jac->as_amg_alpha_opts[3] = 6; 1475 jac->as_amg_alpha_opts[4] = 4; 1476 jac->as_amg_alpha_theta = 0.25; 1477 /* Vector Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */ 1478 jac->as_amg_beta_opts[0] = 10; 1479 jac->as_amg_beta_opts[1] = 1; 1480 jac->as_amg_beta_opts[2] = 6; 1481 jac->as_amg_beta_opts[3] = 6; 1482 jac->as_amg_beta_opts[4] = 4; 1483 jac->as_amg_beta_theta = 0.25; 1484 PetscStackCallStandard(HYPRE_ADSSetPrintLevel,(jac->hsolver,jac->as_print)); 1485 PetscStackCallStandard(HYPRE_ADSSetMaxIter,(jac->hsolver,jac->as_max_iter)); 1486 PetscStackCallStandard(HYPRE_ADSSetCycleType,(jac->hsolver,jac->ams_cycle_type)); 1487 PetscStackCallStandard(HYPRE_ADSSetTol,(jac->hsolver,jac->as_tol)); 1488 PetscStackCallStandard(HYPRE_ADSSetSmoothingOptions,(jac->hsolver,jac->as_relax_type, 1489 jac->as_relax_times, 1490 jac->as_relax_weight, 1491 jac->as_omega)); 1492 PetscStackCallStandard(HYPRE_ADSSetAMSOptions,(jac->hsolver,jac->ams_cycle_type, /* AMG coarsen type */ 1493 jac->as_amg_alpha_opts[0], /* AMG coarsen type */ 1494 jac->as_amg_alpha_opts[1], /* AMG agg_levels */ 1495 jac->as_amg_alpha_opts[2], /* AMG relax_type */ 1496 jac->as_amg_alpha_theta, 1497 jac->as_amg_alpha_opts[3], /* AMG interp_type */ 1498 jac->as_amg_alpha_opts[4])); /* AMG Pmax */ 1499 PetscStackCallStandard(HYPRE_ADSSetAMGOptions,(jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */ 1500 jac->as_amg_beta_opts[1], /* AMG agg_levels */ 1501 jac->as_amg_beta_opts[2], /* AMG relax_type */ 1502 jac->as_amg_beta_theta, 1503 jac->as_amg_beta_opts[3], /* AMG interp_type */ 1504 jac->as_amg_beta_opts[4])); /* AMG Pmax */ 1505 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCSetCoordinates_C",PCSetCoordinates_HYPRE);CHKERRQ(ierr); 1506 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteGradient_C",PCHYPRESetDiscreteGradient_HYPRE);CHKERRQ(ierr); 1507 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteCurl_C",PCHYPRESetDiscreteCurl_HYPRE);CHKERRQ(ierr); 1508 PetscFunctionReturn(0); 1509 } 1510 ierr = PetscFree(jac->hypre_type);CHKERRQ(ierr); 1511 1512 jac->hypre_type = NULL; 1513 SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown HYPRE preconditioner %s; Choices are pilut, parasails, boomeramg, ams",name); 1514 PetscFunctionReturn(0); 1515 } 1516 1517 /* 1518 It only gets here if the HYPRE type has not been set before the call to 1519 ...SetFromOptions() which actually is most of the time 1520 */ 1521 #undef __FUNCT__ 1522 #define __FUNCT__ "PCSetFromOptions_HYPRE" 1523 static PetscErrorCode PCSetFromOptions_HYPRE(PetscOptions *PetscOptionsObject,PC pc) 1524 { 1525 PetscErrorCode ierr; 1526 PetscInt indx; 1527 const char *type[] = {"pilut","parasails","boomeramg","ams","ads"}; 1528 PetscBool flg; 1529 1530 PetscFunctionBegin; 1531 ierr = PetscOptionsHead(PetscOptionsObject,"HYPRE preconditioner options");CHKERRQ(ierr); 1532 ierr = PetscOptionsEList("-pc_hypre_type","HYPRE preconditioner type","PCHYPRESetType",type,4,"boomeramg",&indx,&flg);CHKERRQ(ierr); 1533 if (flg) { 1534 ierr = PCHYPRESetType_HYPRE(pc,type[indx]);CHKERRQ(ierr); 1535 } else { 1536 ierr = PCHYPRESetType_HYPRE(pc,"boomeramg");CHKERRQ(ierr); 1537 } 1538 if (pc->ops->setfromoptions) { 1539 ierr = pc->ops->setfromoptions(PetscOptionsObject,pc);CHKERRQ(ierr); 1540 } 1541 ierr = PetscOptionsTail();CHKERRQ(ierr); 1542 PetscFunctionReturn(0); 1543 } 1544 1545 #undef __FUNCT__ 1546 #define __FUNCT__ "PCHYPRESetType" 1547 /*@C 1548 PCHYPRESetType - Sets which hypre preconditioner you wish to use 1549 1550 Input Parameters: 1551 + pc - the preconditioner context 1552 - name - either pilut, parasails, boomeramg, ams, ads 1553 1554 Options Database Keys: 1555 -pc_hypre_type - One of pilut, parasails, boomeramg, ams, ads 1556 1557 Level: intermediate 1558 1559 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 1560 PCHYPRE 1561 1562 @*/ 1563 PetscErrorCode PCHYPRESetType(PC pc,const char name[]) 1564 { 1565 PetscErrorCode ierr; 1566 1567 PetscFunctionBegin; 1568 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1569 PetscValidCharPointer(name,2); 1570 ierr = PetscTryMethod(pc,"PCHYPRESetType_C",(PC,const char[]),(pc,name));CHKERRQ(ierr); 1571 PetscFunctionReturn(0); 1572 } 1573 1574 #undef __FUNCT__ 1575 #define __FUNCT__ "PCHYPREGetType" 1576 /*@C 1577 PCHYPREGetType - Gets which hypre preconditioner you are using 1578 1579 Input Parameter: 1580 . pc - the preconditioner context 1581 1582 Output Parameter: 1583 . name - either pilut, parasails, boomeramg, ams, ads 1584 1585 Level: intermediate 1586 1587 .seealso: PCCreate(), PCHYPRESetType(), PCType (for list of available types), PC, 1588 PCHYPRE 1589 1590 @*/ 1591 PetscErrorCode PCHYPREGetType(PC pc,const char *name[]) 1592 { 1593 PetscErrorCode ierr; 1594 1595 PetscFunctionBegin; 1596 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1597 PetscValidPointer(name,2); 1598 ierr = PetscTryMethod(pc,"PCHYPREGetType_C",(PC,const char*[]),(pc,name));CHKERRQ(ierr); 1599 PetscFunctionReturn(0); 1600 } 1601 1602 /*MC 1603 PCHYPRE - Allows you to use the matrix element based preconditioners in the LLNL package hypre 1604 1605 Options Database Keys: 1606 + -pc_hypre_type - One of pilut, parasails, boomeramg, ams, ads 1607 - Too many others to list, run with -pc_type hypre -pc_hypre_type XXX -help to see options for the XXX 1608 preconditioner 1609 1610 Level: intermediate 1611 1612 Notes: Apart from pc_hypre_type (for which there is PCHYPRESetType()), 1613 the many hypre options can ONLY be set via the options database (e.g. the command line 1614 or with PetscOptionsSetValue(), there are no functions to set them) 1615 1616 The options -pc_hypre_boomeramg_max_iter and -pc_hypre_boomeramg_rtol refer to the number of iterations 1617 (V-cycles) and tolerance that boomeramg does EACH time it is called. So for example, if 1618 -pc_hypre_boomeramg_max_iter is set to 2 then 2-V-cycles are being used to define the preconditioner 1619 (-pc_hypre_boomeramg_rtol should be set to 0.0 - the default - to strictly use a fixed number of 1620 iterations per hypre call). -ksp_max_it and -ksp_rtol STILL determine the total number of iterations 1621 and tolerance for the Krylov solver. For example, if -pc_hypre_boomeramg_max_iter is 2 and -ksp_max_it is 10 1622 then AT MOST twenty V-cycles of boomeramg will be called. 1623 1624 Note that the option -pc_hypre_boomeramg_relax_type_all defaults to symmetric relaxation 1625 (symmetric-SOR/Jacobi), which is required for Krylov solvers like CG that expect symmetry. 1626 Otherwise, you may want to use -pc_hypre_boomeramg_relax_type_all SOR/Jacobi. 1627 If you wish to use BoomerAMG WITHOUT a Krylov method use -ksp_type richardson NOT -ksp_type preonly 1628 and use -ksp_max_it to control the number of V-cycles. 1629 (see the PETSc FAQ.html at the PETSc website under the Documentation tab). 1630 1631 2007-02-03 Using HYPRE-1.11.1b, the routine HYPRE_BoomerAMGSolveT and the option 1632 -pc_hypre_parasails_reuse were failing with SIGSEGV. Dalcin L. 1633 1634 See PCPFMG for access to the hypre Struct PFMG solver 1635 1636 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 1637 PCHYPRESetType(), PCPFMG 1638 1639 M*/ 1640 1641 #undef __FUNCT__ 1642 #define __FUNCT__ "PCCreate_HYPRE" 1643 PETSC_EXTERN PetscErrorCode PCCreate_HYPRE(PC pc) 1644 { 1645 PC_HYPRE *jac; 1646 PetscErrorCode ierr; 1647 1648 PetscFunctionBegin; 1649 ierr = PetscNewLog(pc,&jac);CHKERRQ(ierr); 1650 1651 pc->data = jac; 1652 pc->ops->destroy = PCDestroy_HYPRE; 1653 pc->ops->setfromoptions = PCSetFromOptions_HYPRE; 1654 pc->ops->setup = PCSetUp_HYPRE; 1655 pc->ops->apply = PCApply_HYPRE; 1656 jac->comm_hypre = MPI_COMM_NULL; 1657 jac->hypre_type = NULL; 1658 jac->coords[0] = NULL; 1659 jac->coords[1] = NULL; 1660 jac->coords[2] = NULL; 1661 jac->constants[0] = NULL; 1662 jac->constants[1] = NULL; 1663 jac->constants[2] = NULL; 1664 jac->G = NULL; 1665 jac->C = NULL; 1666 jac->alpha_Poisson = NULL; 1667 jac->beta_Poisson = NULL; 1668 jac->setdim = NULL; 1669 /* duplicate communicator for hypre */ 1670 ierr = MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(jac->comm_hypre));CHKERRQ(ierr); 1671 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetType_C",PCHYPRESetType_HYPRE);CHKERRQ(ierr); 1672 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCHYPREGetType_C",PCHYPREGetType_HYPRE);CHKERRQ(ierr); 1673 PetscFunctionReturn(0); 1674 } 1675 1676 /* ---------------------------------------------------------------------------------------------------------------------------------*/ 1677 1678 /* this include is needed ONLY to allow access to the private data inside the Mat object specific to hypre */ 1679 #include <petsc/private/matimpl.h> 1680 1681 typedef struct { 1682 MPI_Comm hcomm; /* does not share comm with HYPRE_StructMatrix because need to create solver before getting matrix */ 1683 HYPRE_StructSolver hsolver; 1684 1685 /* keep copy of PFMG options used so may view them */ 1686 PetscInt its; 1687 double tol; 1688 PetscInt relax_type; 1689 PetscInt rap_type; 1690 PetscInt num_pre_relax,num_post_relax; 1691 PetscInt max_levels; 1692 } PC_PFMG; 1693 1694 #undef __FUNCT__ 1695 #define __FUNCT__ "PCDestroy_PFMG" 1696 PetscErrorCode PCDestroy_PFMG(PC pc) 1697 { 1698 PetscErrorCode ierr; 1699 PC_PFMG *ex = (PC_PFMG*) pc->data; 1700 1701 PetscFunctionBegin; 1702 if (ex->hsolver) PetscStackCallStandard(HYPRE_StructPFMGDestroy,(ex->hsolver)); 1703 ierr = MPI_Comm_free(&ex->hcomm);CHKERRQ(ierr); 1704 ierr = PetscFree(pc->data);CHKERRQ(ierr); 1705 PetscFunctionReturn(0); 1706 } 1707 1708 static const char *PFMGRelaxType[] = {"Jacobi","Weighted-Jacobi","symmetric-Red/Black-Gauss-Seidel","Red/Black-Gauss-Seidel"}; 1709 static const char *PFMGRAPType[] = {"Galerkin","non-Galerkin"}; 1710 1711 #undef __FUNCT__ 1712 #define __FUNCT__ "PCView_PFMG" 1713 PetscErrorCode PCView_PFMG(PC pc,PetscViewer viewer) 1714 { 1715 PetscErrorCode ierr; 1716 PetscBool iascii; 1717 PC_PFMG *ex = (PC_PFMG*) pc->data; 1718 1719 PetscFunctionBegin; 1720 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 1721 if (iascii) { 1722 ierr = PetscViewerASCIIPrintf(viewer," HYPRE PFMG preconditioning\n");CHKERRQ(ierr); 1723 ierr = PetscViewerASCIIPrintf(viewer," HYPRE PFMG: max iterations %d\n",ex->its);CHKERRQ(ierr); 1724 ierr = PetscViewerASCIIPrintf(viewer," HYPRE PFMG: tolerance %g\n",ex->tol);CHKERRQ(ierr); 1725 ierr = PetscViewerASCIIPrintf(viewer," HYPRE PFMG: relax type %s\n",PFMGRelaxType[ex->relax_type]);CHKERRQ(ierr); 1726 ierr = PetscViewerASCIIPrintf(viewer," HYPRE PFMG: RAP type %s\n",PFMGRAPType[ex->rap_type]);CHKERRQ(ierr); 1727 ierr = PetscViewerASCIIPrintf(viewer," HYPRE PFMG: number pre-relax %d post-relax %d\n",ex->num_pre_relax,ex->num_post_relax);CHKERRQ(ierr); 1728 ierr = PetscViewerASCIIPrintf(viewer," HYPRE PFMG: max levels %d\n",ex->max_levels);CHKERRQ(ierr); 1729 } 1730 PetscFunctionReturn(0); 1731 } 1732 1733 1734 #undef __FUNCT__ 1735 #define __FUNCT__ "PCSetFromOptions_PFMG" 1736 PetscErrorCode PCSetFromOptions_PFMG(PetscOptions *PetscOptionsObject,PC pc) 1737 { 1738 PetscErrorCode ierr; 1739 PC_PFMG *ex = (PC_PFMG*) pc->data; 1740 PetscBool flg = PETSC_FALSE; 1741 1742 PetscFunctionBegin; 1743 ierr = PetscOptionsHead(PetscOptionsObject,"PFMG options");CHKERRQ(ierr); 1744 ierr = PetscOptionsBool("-pc_pfmg_print_statistics","Print statistics","HYPRE_StructPFMGSetPrintLevel",flg,&flg,NULL);CHKERRQ(ierr); 1745 if (flg) { 1746 int level=3; 1747 PetscStackCallStandard(HYPRE_StructPFMGSetPrintLevel,(ex->hsolver,level)); 1748 } 1749 ierr = PetscOptionsInt("-pc_pfmg_its","Number of iterations of PFMG to use as preconditioner","HYPRE_StructPFMGSetMaxIter",ex->its,&ex->its,NULL);CHKERRQ(ierr); 1750 PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,(ex->hsolver,ex->its)); 1751 ierr = PetscOptionsInt("-pc_pfmg_num_pre_relax","Number of smoothing steps before coarse grid","HYPRE_StructPFMGSetNumPreRelax",ex->num_pre_relax,&ex->num_pre_relax,NULL);CHKERRQ(ierr); 1752 PetscStackCallStandard(HYPRE_StructPFMGSetNumPreRelax,(ex->hsolver,ex->num_pre_relax)); 1753 ierr = PetscOptionsInt("-pc_pfmg_num_post_relax","Number of smoothing steps after coarse grid","HYPRE_StructPFMGSetNumPostRelax",ex->num_post_relax,&ex->num_post_relax,NULL);CHKERRQ(ierr); 1754 PetscStackCallStandard(HYPRE_StructPFMGSetNumPostRelax,(ex->hsolver,ex->num_post_relax)); 1755 1756 ierr = PetscOptionsInt("-pc_pfmg_max_levels","Max Levels for MG hierarchy","HYPRE_StructPFMGSetMaxLevels",ex->max_levels,&ex->max_levels,NULL);CHKERRQ(ierr); 1757 PetscStackCallStandard(HYPRE_StructPFMGSetMaxLevels,(ex->hsolver,ex->max_levels)); 1758 1759 ierr = PetscOptionsReal("-pc_pfmg_tol","Tolerance of PFMG","HYPRE_StructPFMGSetTol",ex->tol,&ex->tol,NULL);CHKERRQ(ierr); 1760 PetscStackCallStandard(HYPRE_StructPFMGSetTol,(ex->hsolver,ex->tol)); 1761 ierr = PetscOptionsEList("-pc_pfmg_relax_type","Relax type for the up and down cycles","HYPRE_StructPFMGSetRelaxType",PFMGRelaxType,ALEN(PFMGRelaxType),PFMGRelaxType[ex->relax_type],&ex->relax_type,NULL);CHKERRQ(ierr); 1762 PetscStackCallStandard(HYPRE_StructPFMGSetRelaxType,(ex->hsolver, ex->relax_type)); 1763 ierr = PetscOptionsEList("-pc_pfmg_rap_type","RAP type","HYPRE_StructPFMGSetRAPType",PFMGRAPType,ALEN(PFMGRAPType),PFMGRAPType[ex->rap_type],&ex->rap_type,NULL);CHKERRQ(ierr); 1764 PetscStackCallStandard(HYPRE_StructPFMGSetRAPType,(ex->hsolver, ex->rap_type)); 1765 ierr = PetscOptionsTail();CHKERRQ(ierr); 1766 PetscFunctionReturn(0); 1767 } 1768 1769 #undef __FUNCT__ 1770 #define __FUNCT__ "PCApply_PFMG" 1771 PetscErrorCode PCApply_PFMG(PC pc,Vec x,Vec y) 1772 { 1773 PetscErrorCode ierr; 1774 PC_PFMG *ex = (PC_PFMG*) pc->data; 1775 PetscScalar *yy; 1776 const PetscScalar *xx; 1777 PetscInt ilower[3],iupper[3]; 1778 Mat_HYPREStruct *mx = (Mat_HYPREStruct*)(pc->pmat->data); 1779 1780 PetscFunctionBegin; 1781 ierr = PetscCitationsRegister(hypreCitation,&cite);CHKERRQ(ierr); 1782 ierr = DMDAGetCorners(mx->da,&ilower[0],&ilower[1],&ilower[2],&iupper[0],&iupper[1],&iupper[2]);CHKERRQ(ierr); 1783 iupper[0] += ilower[0] - 1; 1784 iupper[1] += ilower[1] - 1; 1785 iupper[2] += ilower[2] - 1; 1786 1787 /* copy x values over to hypre */ 1788 PetscStackCallStandard(HYPRE_StructVectorSetConstantValues,(mx->hb,0.0)); 1789 ierr = VecGetArrayRead(x,&xx);CHKERRQ(ierr); 1790 PetscStackCallStandard(HYPRE_StructVectorSetBoxValues,(mx->hb,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,(PetscScalar*)xx)); 1791 ierr = VecRestoreArrayRead(x,&xx);CHKERRQ(ierr); 1792 PetscStackCallStandard(HYPRE_StructVectorAssemble,(mx->hb)); 1793 PetscStackCallStandard(HYPRE_StructPFMGSolve,(ex->hsolver,mx->hmat,mx->hb,mx->hx)); 1794 1795 /* copy solution values back to PETSc */ 1796 ierr = VecGetArray(y,&yy);CHKERRQ(ierr); 1797 PetscStackCallStandard(HYPRE_StructVectorGetBoxValues,(mx->hx,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,yy)); 1798 ierr = VecRestoreArray(y,&yy);CHKERRQ(ierr); 1799 PetscFunctionReturn(0); 1800 } 1801 1802 #undef __FUNCT__ 1803 #define __FUNCT__ "PCApplyRichardson_PFMG" 1804 static PetscErrorCode PCApplyRichardson_PFMG(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscBool guesszero,PetscInt *outits,PCRichardsonConvergedReason *reason) 1805 { 1806 PC_PFMG *jac = (PC_PFMG*)pc->data; 1807 PetscErrorCode ierr; 1808 PetscInt oits; 1809 1810 PetscFunctionBegin; 1811 ierr = PetscCitationsRegister(hypreCitation,&cite);CHKERRQ(ierr); 1812 PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,(jac->hsolver,its*jac->its)); 1813 PetscStackCallStandard(HYPRE_StructPFMGSetTol,(jac->hsolver,rtol)); 1814 1815 ierr = PCApply_PFMG(pc,b,y);CHKERRQ(ierr); 1816 PetscStackCallStandard(HYPRE_StructPFMGGetNumIterations,(jac->hsolver,(HYPRE_Int *)&oits)); 1817 *outits = oits; 1818 if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS; 1819 else *reason = PCRICHARDSON_CONVERGED_RTOL; 1820 PetscStackCallStandard(HYPRE_StructPFMGSetTol,(jac->hsolver,jac->tol)); 1821 PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,(jac->hsolver,jac->its)); 1822 PetscFunctionReturn(0); 1823 } 1824 1825 1826 #undef __FUNCT__ 1827 #define __FUNCT__ "PCSetUp_PFMG" 1828 PetscErrorCode PCSetUp_PFMG(PC pc) 1829 { 1830 PetscErrorCode ierr; 1831 PC_PFMG *ex = (PC_PFMG*) pc->data; 1832 Mat_HYPREStruct *mx = (Mat_HYPREStruct*)(pc->pmat->data); 1833 PetscBool flg; 1834 1835 PetscFunctionBegin; 1836 ierr = PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRESTRUCT,&flg);CHKERRQ(ierr); 1837 if (!flg) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Must use MATHYPRESTRUCT with this preconditioner"); 1838 1839 /* create the hypre solver object and set its information */ 1840 if (ex->hsolver) PetscStackCallStandard(HYPRE_StructPFMGDestroy,(ex->hsolver)); 1841 PetscStackCallStandard(HYPRE_StructPFMGCreate,(ex->hcomm,&ex->hsolver)); 1842 PetscStackCallStandard(HYPRE_StructPFMGSetup,(ex->hsolver,mx->hmat,mx->hb,mx->hx)); 1843 PetscStackCallStandard(HYPRE_StructPFMGSetZeroGuess,(ex->hsolver)); 1844 PetscFunctionReturn(0); 1845 } 1846 1847 1848 /*MC 1849 PCPFMG - the hypre PFMG multigrid solver 1850 1851 Level: advanced 1852 1853 Options Database: 1854 + -pc_pfmg_its <its> number of iterations of PFMG to use as preconditioner 1855 . -pc_pfmg_num_pre_relax <steps> number of smoothing steps before coarse grid 1856 . -pc_pfmg_num_post_relax <steps> number of smoothing steps after coarse grid 1857 . -pc_pfmg_tol <tol> tolerance of PFMG 1858 . -pc_pfmg_relax_type -relaxation type for the up and down cycles, one of Jacobi,Weighted-Jacobi,symmetric-Red/Black-Gauss-Seidel,Red/Black-Gauss-Seidel 1859 - -pc_pfmg_rap_type - type of coarse matrix generation, one of Galerkin,non-Galerkin 1860 1861 Notes: This is for CELL-centered descretizations 1862 1863 This must be used with the MATHYPRESTRUCT matrix type. 1864 This is less general than in hypre, it supports only one block per process defined by a PETSc DMDA. 1865 1866 .seealso: PCMG, MATHYPRESTRUCT 1867 M*/ 1868 1869 #undef __FUNCT__ 1870 #define __FUNCT__ "PCCreate_PFMG" 1871 PETSC_EXTERN PetscErrorCode PCCreate_PFMG(PC pc) 1872 { 1873 PetscErrorCode ierr; 1874 PC_PFMG *ex; 1875 1876 PetscFunctionBegin; 1877 ierr = PetscNew(&ex);CHKERRQ(ierr); \ 1878 pc->data = ex; 1879 1880 ex->its = 1; 1881 ex->tol = 1.e-8; 1882 ex->relax_type = 1; 1883 ex->rap_type = 0; 1884 ex->num_pre_relax = 1; 1885 ex->num_post_relax = 1; 1886 ex->max_levels = 0; 1887 1888 pc->ops->setfromoptions = PCSetFromOptions_PFMG; 1889 pc->ops->view = PCView_PFMG; 1890 pc->ops->destroy = PCDestroy_PFMG; 1891 pc->ops->apply = PCApply_PFMG; 1892 pc->ops->applyrichardson = PCApplyRichardson_PFMG; 1893 pc->ops->setup = PCSetUp_PFMG; 1894 1895 ierr = MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(ex->hcomm));CHKERRQ(ierr); 1896 PetscStackCallStandard(HYPRE_StructPFMGCreate,(ex->hcomm,&ex->hsolver)); 1897 PetscFunctionReturn(0); 1898 } 1899 1900 /* ---------------------------------------------------------------------------------------------------------------------------------------------------*/ 1901 1902 /* we know we are working with a HYPRE_SStructMatrix */ 1903 typedef struct { 1904 MPI_Comm hcomm; /* does not share comm with HYPRE_SStructMatrix because need to create solver before getting matrix */ 1905 HYPRE_SStructSolver ss_solver; 1906 1907 /* keep copy of SYSPFMG options used so may view them */ 1908 PetscInt its; 1909 double tol; 1910 PetscInt relax_type; 1911 PetscInt num_pre_relax,num_post_relax; 1912 } PC_SysPFMG; 1913 1914 #undef __FUNCT__ 1915 #define __FUNCT__ "PCDestroy_SysPFMG" 1916 PetscErrorCode PCDestroy_SysPFMG(PC pc) 1917 { 1918 PetscErrorCode ierr; 1919 PC_SysPFMG *ex = (PC_SysPFMG*) pc->data; 1920 1921 PetscFunctionBegin; 1922 if (ex->ss_solver) PetscStackCallStandard(HYPRE_SStructSysPFMGDestroy,(ex->ss_solver)); 1923 ierr = MPI_Comm_free(&ex->hcomm);CHKERRQ(ierr); 1924 ierr = PetscFree(pc->data);CHKERRQ(ierr); 1925 PetscFunctionReturn(0); 1926 } 1927 1928 static const char *SysPFMGRelaxType[] = {"Weighted-Jacobi","Red/Black-Gauss-Seidel"}; 1929 1930 #undef __FUNCT__ 1931 #define __FUNCT__ "PCView_SysPFMG" 1932 PetscErrorCode PCView_SysPFMG(PC pc,PetscViewer viewer) 1933 { 1934 PetscErrorCode ierr; 1935 PetscBool iascii; 1936 PC_SysPFMG *ex = (PC_SysPFMG*) pc->data; 1937 1938 PetscFunctionBegin; 1939 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 1940 if (iascii) { 1941 ierr = PetscViewerASCIIPrintf(viewer," HYPRE SysPFMG preconditioning\n");CHKERRQ(ierr); 1942 ierr = PetscViewerASCIIPrintf(viewer," HYPRE SysPFMG: max iterations %d\n",ex->its);CHKERRQ(ierr); 1943 ierr = PetscViewerASCIIPrintf(viewer," HYPRE SysPFMG: tolerance %g\n",ex->tol);CHKERRQ(ierr); 1944 ierr = PetscViewerASCIIPrintf(viewer," HYPRE SysPFMG: relax type %s\n",PFMGRelaxType[ex->relax_type]);CHKERRQ(ierr); 1945 ierr = PetscViewerASCIIPrintf(viewer," HYPRE SysPFMG: number pre-relax %d post-relax %d\n",ex->num_pre_relax,ex->num_post_relax);CHKERRQ(ierr); 1946 } 1947 PetscFunctionReturn(0); 1948 } 1949 1950 1951 #undef __FUNCT__ 1952 #define __FUNCT__ "PCSetFromOptions_SysPFMG" 1953 PetscErrorCode PCSetFromOptions_SysPFMG(PetscOptions *PetscOptionsObject,PC pc) 1954 { 1955 PetscErrorCode ierr; 1956 PC_SysPFMG *ex = (PC_SysPFMG*) pc->data; 1957 PetscBool flg = PETSC_FALSE; 1958 1959 PetscFunctionBegin; 1960 ierr = PetscOptionsHead(PetscOptionsObject,"SysPFMG options");CHKERRQ(ierr); 1961 ierr = PetscOptionsBool("-pc_syspfmg_print_statistics","Print statistics","HYPRE_SStructSysPFMGSetPrintLevel",flg,&flg,NULL);CHKERRQ(ierr); 1962 if (flg) { 1963 int level=3; 1964 PetscStackCallStandard(HYPRE_SStructSysPFMGSetPrintLevel,(ex->ss_solver,level)); 1965 } 1966 ierr = PetscOptionsInt("-pc_syspfmg_its","Number of iterations of SysPFMG to use as preconditioner","HYPRE_SStructSysPFMGSetMaxIter",ex->its,&ex->its,NULL);CHKERRQ(ierr); 1967 PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,(ex->ss_solver,ex->its)); 1968 ierr = PetscOptionsInt("-pc_syspfmg_num_pre_relax","Number of smoothing steps before coarse grid","HYPRE_SStructSysPFMGSetNumPreRelax",ex->num_pre_relax,&ex->num_pre_relax,NULL);CHKERRQ(ierr); 1969 PetscStackCallStandard(HYPRE_SStructSysPFMGSetNumPreRelax,(ex->ss_solver,ex->num_pre_relax)); 1970 ierr = PetscOptionsInt("-pc_syspfmg_num_post_relax","Number of smoothing steps after coarse grid","HYPRE_SStructSysPFMGSetNumPostRelax",ex->num_post_relax,&ex->num_post_relax,NULL);CHKERRQ(ierr); 1971 PetscStackCallStandard(HYPRE_SStructSysPFMGSetNumPostRelax,(ex->ss_solver,ex->num_post_relax)); 1972 1973 ierr = PetscOptionsReal("-pc_syspfmg_tol","Tolerance of SysPFMG","HYPRE_SStructSysPFMGSetTol",ex->tol,&ex->tol,NULL);CHKERRQ(ierr); 1974 PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,(ex->ss_solver,ex->tol)); 1975 ierr = PetscOptionsEList("-pc_syspfmg_relax_type","Relax type for the up and down cycles","HYPRE_SStructSysPFMGSetRelaxType",SysPFMGRelaxType,4,SysPFMGRelaxType[ex->relax_type],&ex->relax_type,NULL);CHKERRQ(ierr); 1976 PetscStackCallStandard(HYPRE_SStructSysPFMGSetRelaxType,(ex->ss_solver, ex->relax_type)); 1977 ierr = PetscOptionsTail();CHKERRQ(ierr); 1978 PetscFunctionReturn(0); 1979 } 1980 1981 #undef __FUNCT__ 1982 #define __FUNCT__ "PCApply_SysPFMG" 1983 PetscErrorCode PCApply_SysPFMG(PC pc,Vec x,Vec y) 1984 { 1985 PetscErrorCode ierr; 1986 PC_SysPFMG *ex = (PC_SysPFMG*) pc->data; 1987 PetscScalar *yy; 1988 const PetscScalar *xx; 1989 PetscInt ilower[3],iupper[3]; 1990 Mat_HYPRESStruct *mx = (Mat_HYPRESStruct*)(pc->pmat->data); 1991 PetscInt ordering= mx->dofs_order; 1992 PetscInt nvars = mx->nvars; 1993 PetscInt part = 0; 1994 PetscInt size; 1995 PetscInt i; 1996 1997 PetscFunctionBegin; 1998 ierr = PetscCitationsRegister(hypreCitation,&cite);CHKERRQ(ierr); 1999 ierr = DMDAGetCorners(mx->da,&ilower[0],&ilower[1],&ilower[2],&iupper[0],&iupper[1],&iupper[2]);CHKERRQ(ierr); 2000 iupper[0] += ilower[0] - 1; 2001 iupper[1] += ilower[1] - 1; 2002 iupper[2] += ilower[2] - 1; 2003 2004 size = 1; 2005 for (i= 0; i< 3; i++) size *= (iupper[i]-ilower[i]+1); 2006 2007 /* copy x values over to hypre for variable ordering */ 2008 if (ordering) { 2009 PetscStackCallStandard(HYPRE_SStructVectorSetConstantValues,(mx->ss_b,0.0)); 2010 ierr = VecGetArrayRead(x,&xx);CHKERRQ(ierr); 2011 for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorSetBoxValues,(mx->ss_b,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,(PetscScalar*)xx+(size*i))); 2012 ierr = VecRestoreArrayRead(x,&xx);CHKERRQ(ierr); 2013 PetscStackCallStandard(HYPRE_SStructVectorAssemble,(mx->ss_b)); 2014 PetscStackCallStandard(HYPRE_SStructMatrixMatvec,(1.0,mx->ss_mat,mx->ss_b,0.0,mx->ss_x)); 2015 PetscStackCallStandard(HYPRE_SStructSysPFMGSolve,(ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x)); 2016 2017 /* copy solution values back to PETSc */ 2018 ierr = VecGetArray(y,&yy);CHKERRQ(ierr); 2019 for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorGetBoxValues,(mx->ss_x,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,yy+(size*i))); 2020 ierr = VecRestoreArray(y,&yy);CHKERRQ(ierr); 2021 } else { /* nodal ordering must be mapped to variable ordering for sys_pfmg */ 2022 PetscScalar *z; 2023 PetscInt j, k; 2024 2025 ierr = PetscMalloc1(nvars*size,&z);CHKERRQ(ierr); 2026 PetscStackCallStandard(HYPRE_SStructVectorSetConstantValues,(mx->ss_b,0.0)); 2027 ierr = VecGetArrayRead(x,&xx);CHKERRQ(ierr); 2028 2029 /* transform nodal to hypre's variable ordering for sys_pfmg */ 2030 for (i= 0; i< size; i++) { 2031 k= i*nvars; 2032 for (j= 0; j< nvars; j++) z[j*size+i]= xx[k+j]; 2033 } 2034 for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorSetBoxValues,(mx->ss_b,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,z+(size*i))); 2035 ierr = VecRestoreArrayRead(x,&xx);CHKERRQ(ierr); 2036 PetscStackCallStandard(HYPRE_SStructVectorAssemble,(mx->ss_b)); 2037 PetscStackCallStandard(HYPRE_SStructSysPFMGSolve,(ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x)); 2038 2039 /* copy solution values back to PETSc */ 2040 ierr = VecGetArray(y,&yy);CHKERRQ(ierr); 2041 for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorGetBoxValues,(mx->ss_x,part,(HYPRE_Int *)ilower,(HYPRE_Int *)iupper,i,z+(size*i))); 2042 /* transform hypre's variable ordering for sys_pfmg to nodal ordering */ 2043 for (i= 0; i< size; i++) { 2044 k= i*nvars; 2045 for (j= 0; j< nvars; j++) yy[k+j]= z[j*size+i]; 2046 } 2047 ierr = VecRestoreArray(y,&yy);CHKERRQ(ierr); 2048 ierr = PetscFree(z);CHKERRQ(ierr); 2049 } 2050 PetscFunctionReturn(0); 2051 } 2052 2053 #undef __FUNCT__ 2054 #define __FUNCT__ "PCApplyRichardson_SysPFMG" 2055 static PetscErrorCode PCApplyRichardson_SysPFMG(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscBool guesszero,PetscInt *outits,PCRichardsonConvergedReason *reason) 2056 { 2057 PC_SysPFMG *jac = (PC_SysPFMG*)pc->data; 2058 PetscErrorCode ierr; 2059 PetscInt oits; 2060 2061 PetscFunctionBegin; 2062 ierr = PetscCitationsRegister(hypreCitation,&cite);CHKERRQ(ierr); 2063 PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,(jac->ss_solver,its*jac->its)); 2064 PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,(jac->ss_solver,rtol)); 2065 ierr = PCApply_SysPFMG(pc,b,y);CHKERRQ(ierr); 2066 PetscStackCallStandard(HYPRE_SStructSysPFMGGetNumIterations,(jac->ss_solver,(HYPRE_Int *)&oits)); 2067 *outits = oits; 2068 if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS; 2069 else *reason = PCRICHARDSON_CONVERGED_RTOL; 2070 PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,(jac->ss_solver,jac->tol)); 2071 PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,(jac->ss_solver,jac->its)); 2072 PetscFunctionReturn(0); 2073 } 2074 2075 2076 #undef __FUNCT__ 2077 #define __FUNCT__ "PCSetUp_SysPFMG" 2078 PetscErrorCode PCSetUp_SysPFMG(PC pc) 2079 { 2080 PetscErrorCode ierr; 2081 PC_SysPFMG *ex = (PC_SysPFMG*) pc->data; 2082 Mat_HYPRESStruct *mx = (Mat_HYPRESStruct*)(pc->pmat->data); 2083 PetscBool flg; 2084 2085 PetscFunctionBegin; 2086 ierr = PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRESSTRUCT,&flg);CHKERRQ(ierr); 2087 if (!flg) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Must use MATHYPRESSTRUCT with this preconditioner"); 2088 2089 /* create the hypre sstruct solver object and set its information */ 2090 if (ex->ss_solver) PetscStackCallStandard(HYPRE_SStructSysPFMGDestroy,(ex->ss_solver)); 2091 PetscStackCallStandard(HYPRE_SStructSysPFMGCreate,(ex->hcomm,&ex->ss_solver)); 2092 PetscStackCallStandard(HYPRE_SStructSysPFMGSetZeroGuess,(ex->ss_solver)); 2093 PetscStackCallStandard(HYPRE_SStructSysPFMGSetup,(ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x)); 2094 PetscFunctionReturn(0); 2095 } 2096 2097 2098 /*MC 2099 PCSysPFMG - the hypre SysPFMG multigrid solver 2100 2101 Level: advanced 2102 2103 Options Database: 2104 + -pc_syspfmg_its <its> number of iterations of SysPFMG to use as preconditioner 2105 . -pc_syspfmg_num_pre_relax <steps> number of smoothing steps before coarse grid 2106 . -pc_syspfmg_num_post_relax <steps> number of smoothing steps after coarse grid 2107 . -pc_syspfmg_tol <tol> tolerance of SysPFMG 2108 . -pc_syspfmg_relax_type -relaxation type for the up and down cycles, one of Weighted-Jacobi,Red/Black-Gauss-Seidel 2109 2110 Notes: This is for CELL-centered descretizations 2111 2112 This must be used with the MATHYPRESSTRUCT matrix type. 2113 This is less general than in hypre, it supports only one part, and one block per process defined by a PETSc DMDA. 2114 Also, only cell-centered variables. 2115 2116 .seealso: PCMG, MATHYPRESSTRUCT 2117 M*/ 2118 2119 #undef __FUNCT__ 2120 #define __FUNCT__ "PCCreate_SysPFMG" 2121 PETSC_EXTERN PetscErrorCode PCCreate_SysPFMG(PC pc) 2122 { 2123 PetscErrorCode ierr; 2124 PC_SysPFMG *ex; 2125 2126 PetscFunctionBegin; 2127 ierr = PetscNew(&ex);CHKERRQ(ierr); \ 2128 pc->data = ex; 2129 2130 ex->its = 1; 2131 ex->tol = 1.e-8; 2132 ex->relax_type = 1; 2133 ex->num_pre_relax = 1; 2134 ex->num_post_relax = 1; 2135 2136 pc->ops->setfromoptions = PCSetFromOptions_SysPFMG; 2137 pc->ops->view = PCView_SysPFMG; 2138 pc->ops->destroy = PCDestroy_SysPFMG; 2139 pc->ops->apply = PCApply_SysPFMG; 2140 pc->ops->applyrichardson = PCApplyRichardson_SysPFMG; 2141 pc->ops->setup = PCSetUp_SysPFMG; 2142 2143 ierr = MPI_Comm_dup(PetscObjectComm((PetscObject)pc),&(ex->hcomm));CHKERRQ(ierr); 2144 PetscStackCallStandard(HYPRE_SStructSysPFMGCreate,(ex->hcomm,&ex->ss_solver)); 2145 PetscFunctionReturn(0); 2146 } 2147