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