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