1 /* 2 Provides an interface to the LLNL package hypre 3 */ 4 5 #include <petscpkg_version.h> 6 #include <petsc/private/pcimpl.h> /*I "petscpc.h" I*/ 7 /* this include is needed ONLY to allow access to the private data inside the Mat object specific to hypre */ 8 #include <petsc/private/matimpl.h> 9 #include <petsc/private/vecimpl.h> 10 #include <../src/vec/vec/impls/hypre/vhyp.h> 11 #include <../src/mat/impls/hypre/mhypre.h> 12 #include <../src/dm/impls/da/hypre/mhyp.h> 13 #include <_hypre_parcsr_ls.h> 14 #include <petscmathypre.h> 15 16 #if defined(PETSC_HAVE_HYPRE_DEVICE) 17 #include <petsc/private/deviceimpl.h> 18 #endif 19 20 static PetscBool cite = PETSC_FALSE; 21 static const char hypreCitation[] = "@manual{hypre-web-page,\n title = {{\\sl hypre}: High Performance Preconditioners},\n organization = {Lawrence Livermore National Laboratory},\n note = {\\url{https://computation.llnl.gov/projects/hypre-scalable-linear-solvers-multigrid-methods}}\n}\n"; 22 23 /* 24 Private context (data structure) for the preconditioner. 25 */ 26 typedef struct { 27 HYPRE_Solver hsolver; 28 Mat hpmat; /* MatHYPRE */ 29 30 HYPRE_Int (*destroy)(HYPRE_Solver); 31 HYPRE_Int (*solve)(HYPRE_Solver,HYPRE_ParCSRMatrix,HYPRE_ParVector,HYPRE_ParVector); 32 HYPRE_Int (*setup)(HYPRE_Solver,HYPRE_ParCSRMatrix,HYPRE_ParVector,HYPRE_ParVector); 33 34 MPI_Comm comm_hypre; 35 char *hypre_type; 36 37 /* options for Pilut and BoomerAMG*/ 38 PetscInt maxiter; 39 PetscReal tol; 40 41 /* options for Pilut */ 42 PetscInt factorrowsize; 43 44 /* options for ParaSails */ 45 PetscInt nlevels; 46 PetscReal threshold; 47 PetscReal filter; 48 PetscReal loadbal; 49 PetscInt logging; 50 PetscInt ruse; 51 PetscInt symt; 52 53 /* options for BoomerAMG */ 54 PetscBool printstatistics; 55 56 /* options for BoomerAMG */ 57 PetscInt cycletype; 58 PetscInt maxlevels; 59 PetscReal strongthreshold; 60 PetscReal maxrowsum; 61 PetscInt gridsweeps[3]; 62 PetscInt coarsentype; 63 PetscInt measuretype; 64 PetscInt smoothtype; 65 PetscInt smoothnumlevels; 66 PetscInt eu_level; /* Number of levels for ILU(k) in Euclid */ 67 PetscReal eu_droptolerance; /* Drop tolerance for ILU(k) in Euclid */ 68 PetscInt eu_bj; /* Defines use of Block Jacobi ILU in Euclid */ 69 PetscInt relaxtype[3]; 70 PetscReal relaxweight; 71 PetscReal outerrelaxweight; 72 PetscInt relaxorder; 73 PetscReal truncfactor; 74 PetscBool applyrichardson; 75 PetscInt pmax; 76 PetscInt interptype; 77 PetscInt maxc; 78 PetscInt minc; 79 #if PETSC_PKG_HYPRE_VERSION_GE(2,23,0) 80 char *spgemm_type; // this is a global hypre parameter but is closely associated with BoomerAMG 81 #endif 82 /* GPU */ 83 PetscBool keeptranspose; 84 PetscInt rap2; 85 PetscInt mod_rap2; 86 87 /* AIR */ 88 PetscInt Rtype; 89 PetscReal Rstrongthreshold; 90 PetscReal Rfilterthreshold; 91 PetscInt Adroptype; 92 PetscReal Adroptol; 93 94 PetscInt agg_nl; 95 PetscInt agg_interptype; 96 PetscInt agg_num_paths; 97 PetscBool nodal_relax; 98 PetscInt nodal_relax_levels; 99 100 PetscInt nodal_coarsening; 101 PetscInt nodal_coarsening_diag; 102 PetscInt vec_interp_variant; 103 PetscInt vec_interp_qmax; 104 PetscBool vec_interp_smooth; 105 PetscInt interp_refine; 106 107 /* NearNullSpace support */ 108 VecHYPRE_IJVector *hmnull; 109 HYPRE_ParVector *phmnull; 110 PetscInt n_hmnull; 111 Vec hmnull_constant; 112 113 /* options for AS (Auxiliary Space preconditioners) */ 114 PetscInt as_print; 115 PetscInt as_max_iter; 116 PetscReal as_tol; 117 PetscInt as_relax_type; 118 PetscInt as_relax_times; 119 PetscReal as_relax_weight; 120 PetscReal as_omega; 121 PetscInt as_amg_alpha_opts[5]; /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for vector Poisson (AMS) or Curl problem (ADS) */ 122 PetscReal as_amg_alpha_theta; /* AMG strength for vector Poisson (AMS) or Curl problem (ADS) */ 123 PetscInt as_amg_beta_opts[5]; /* AMG coarsen type, agg_levels, relax_type, interp_type, Pmax for scalar Poisson (AMS) or vector Poisson (ADS) */ 124 PetscReal as_amg_beta_theta; /* AMG strength for scalar Poisson (AMS) or vector Poisson (ADS) */ 125 PetscInt ams_cycle_type; 126 PetscInt ads_cycle_type; 127 128 /* additional data */ 129 Mat G; /* MatHYPRE */ 130 Mat C; /* MatHYPRE */ 131 Mat alpha_Poisson; /* MatHYPRE */ 132 Mat beta_Poisson; /* MatHYPRE */ 133 134 /* extra information for AMS */ 135 PetscInt dim; /* geometrical dimension */ 136 VecHYPRE_IJVector coords[3]; 137 VecHYPRE_IJVector constants[3]; 138 Mat RT_PiFull, RT_Pi[3]; 139 Mat ND_PiFull, ND_Pi[3]; 140 PetscBool ams_beta_is_zero; 141 PetscBool ams_beta_is_zero_part; 142 PetscInt ams_proj_freq; 143 } PC_HYPRE; 144 145 PetscErrorCode PCHYPREGetSolver(PC pc,HYPRE_Solver *hsolver) 146 { 147 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 148 149 PetscFunctionBegin; 150 *hsolver = jac->hsolver; 151 PetscFunctionReturn(0); 152 } 153 154 /* 155 Matrices with AIJ format are created IN PLACE with using (I,J,data) from BoomerAMG. Since the data format in hypre_ParCSRMatrix 156 is different from that used in PETSc, the original hypre_ParCSRMatrix can not be used any more after call this routine. 157 It is used in PCHMG. Other users should avoid using this function. 158 */ 159 static PetscErrorCode PCGetCoarseOperators_BoomerAMG(PC pc,PetscInt *nlevels,Mat *operators[]) 160 { 161 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 162 PetscBool same = PETSC_FALSE; 163 PetscInt num_levels,l; 164 Mat *mattmp; 165 hypre_ParCSRMatrix **A_array; 166 167 PetscFunctionBegin; 168 CHKERRQ(PetscStrcmp(jac->hypre_type,"boomeramg",&same)); 169 PetscCheck(same,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_NOTSAMETYPE,"Hypre type is not BoomerAMG "); 170 num_levels = hypre_ParAMGDataNumLevels((hypre_ParAMGData*) (jac->hsolver)); 171 CHKERRQ(PetscMalloc1(num_levels,&mattmp)); 172 A_array = hypre_ParAMGDataAArray((hypre_ParAMGData*) (jac->hsolver)); 173 for (l=1; l<num_levels; l++) { 174 CHKERRQ(MatCreateFromParCSR(A_array[l],MATAIJ,PETSC_OWN_POINTER, &(mattmp[num_levels-1-l]))); 175 /* We want to own the data, and HYPRE can not touch this matrix any more */ 176 A_array[l] = NULL; 177 } 178 *nlevels = num_levels; 179 *operators = mattmp; 180 PetscFunctionReturn(0); 181 } 182 183 /* 184 Matrices with AIJ format are created IN PLACE with using (I,J,data) from BoomerAMG. Since the data format in hypre_ParCSRMatrix 185 is different from that used in PETSc, the original hypre_ParCSRMatrix can not be used any more after call this routine. 186 It is used in PCHMG. Other users should avoid using this function. 187 */ 188 static PetscErrorCode PCGetInterpolations_BoomerAMG(PC pc,PetscInt *nlevels,Mat *interpolations[]) 189 { 190 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 191 PetscBool same = PETSC_FALSE; 192 PetscInt num_levels,l; 193 Mat *mattmp; 194 hypre_ParCSRMatrix **P_array; 195 196 PetscFunctionBegin; 197 CHKERRQ(PetscStrcmp(jac->hypre_type,"boomeramg",&same)); 198 PetscCheck(same,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_NOTSAMETYPE,"Hypre type is not BoomerAMG "); 199 num_levels = hypre_ParAMGDataNumLevels((hypre_ParAMGData*) (jac->hsolver)); 200 CHKERRQ(PetscMalloc1(num_levels,&mattmp)); 201 P_array = hypre_ParAMGDataPArray((hypre_ParAMGData*) (jac->hsolver)); 202 for (l=1; l<num_levels; l++) { 203 CHKERRQ(MatCreateFromParCSR(P_array[num_levels-1-l],MATAIJ,PETSC_OWN_POINTER, &(mattmp[l-1]))); 204 /* We want to own the data, and HYPRE can not touch this matrix any more */ 205 P_array[num_levels-1-l] = NULL; 206 } 207 *nlevels = num_levels; 208 *interpolations = mattmp; 209 PetscFunctionReturn(0); 210 } 211 212 /* Resets (frees) Hypre's representation of the near null space */ 213 static PetscErrorCode PCHYPREResetNearNullSpace_Private(PC pc) 214 { 215 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 216 PetscInt i; 217 218 PetscFunctionBegin; 219 for (i=0; i<jac->n_hmnull; i++) { 220 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->hmnull[i])); 221 } 222 CHKERRQ(PetscFree(jac->hmnull)); 223 CHKERRQ(PetscFree(jac->phmnull)); 224 CHKERRQ(VecDestroy(&jac->hmnull_constant)); 225 jac->n_hmnull = 0; 226 PetscFunctionReturn(0); 227 } 228 229 static PetscErrorCode PCSetUp_HYPRE(PC pc) 230 { 231 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 232 Mat_HYPRE *hjac; 233 HYPRE_ParCSRMatrix hmat; 234 HYPRE_ParVector bv,xv; 235 PetscBool ishypre; 236 237 PetscFunctionBegin; 238 if (!jac->hypre_type) { 239 CHKERRQ(PCHYPRESetType(pc,"boomeramg")); 240 } 241 242 CHKERRQ(PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRE,&ishypre)); 243 if (!ishypre) { 244 CHKERRQ(MatDestroy(&jac->hpmat)); 245 CHKERRQ(MatConvert(pc->pmat,MATHYPRE,MAT_INITIAL_MATRIX,&jac->hpmat)); 246 CHKERRQ(PetscLogObjectParent((PetscObject)pc,(PetscObject)jac->hpmat)); 247 } else { 248 CHKERRQ(PetscObjectReference((PetscObject)pc->pmat)); 249 CHKERRQ(MatDestroy(&jac->hpmat)); 250 jac->hpmat = pc->pmat; 251 } 252 /* allow debug */ 253 CHKERRQ(MatViewFromOptions(jac->hpmat,NULL,"-pc_hypre_mat_view")); 254 hjac = (Mat_HYPRE*)(jac->hpmat->data); 255 256 /* special case for BoomerAMG */ 257 if (jac->setup == HYPRE_BoomerAMGSetup) { 258 MatNullSpace mnull; 259 PetscBool has_const; 260 PetscInt bs,nvec,i; 261 const Vec *vecs; 262 263 CHKERRQ(MatGetBlockSize(pc->pmat,&bs)); 264 if (bs > 1) PetscStackCallStandard(HYPRE_BoomerAMGSetNumFunctions,jac->hsolver,bs); 265 CHKERRQ(MatGetNearNullSpace(pc->mat, &mnull)); 266 if (mnull) { 267 CHKERRQ(PCHYPREResetNearNullSpace_Private(pc)); 268 CHKERRQ(MatNullSpaceGetVecs(mnull, &has_const, &nvec, &vecs)); 269 CHKERRQ(PetscMalloc1(nvec+1,&jac->hmnull)); 270 CHKERRQ(PetscMalloc1(nvec+1,&jac->phmnull)); 271 for (i=0; i<nvec; i++) { 272 CHKERRQ(VecHYPRE_IJVectorCreate(vecs[i]->map,&jac->hmnull[i])); 273 CHKERRQ(VecHYPRE_IJVectorCopy(vecs[i],jac->hmnull[i])); 274 PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->hmnull[i]->ij,(void**)&jac->phmnull[i]); 275 } 276 if (has_const) { 277 CHKERRQ(MatCreateVecs(pc->pmat,&jac->hmnull_constant,NULL)); 278 CHKERRQ(PetscLogObjectParent((PetscObject)pc,(PetscObject)jac->hmnull_constant)); 279 CHKERRQ(VecSet(jac->hmnull_constant,1)); 280 CHKERRQ(VecNormalize(jac->hmnull_constant,NULL)); 281 CHKERRQ(VecHYPRE_IJVectorCreate(jac->hmnull_constant->map,&jac->hmnull[nvec])); 282 CHKERRQ(VecHYPRE_IJVectorCopy(jac->hmnull_constant,jac->hmnull[nvec])); 283 PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->hmnull[nvec]->ij,(void**)&jac->phmnull[nvec]); 284 nvec++; 285 } 286 PetscStackCallStandard(HYPRE_BoomerAMGSetInterpVectors,jac->hsolver,nvec,jac->phmnull); 287 jac->n_hmnull = nvec; 288 } 289 } 290 291 /* special case for AMS */ 292 if (jac->setup == HYPRE_AMSSetup) { 293 Mat_HYPRE *hm; 294 HYPRE_ParCSRMatrix parcsr; 295 if (!jac->coords[0] && !jac->constants[0] && !(jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1]))) { 296 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE AMS preconditioner needs either the coordinate vectors via PCSetCoordinates() or the edge constant vectors via PCHYPRESetEdgeConstantVectors() or the interpolation matrix via PCHYPRESetInterpolations"); 297 } 298 if (jac->dim) { 299 PetscStackCallStandard(HYPRE_AMSSetDimension,jac->hsolver,jac->dim); 300 } 301 if (jac->constants[0]) { 302 HYPRE_ParVector ozz,zoz,zzo = NULL; 303 PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->constants[0]->ij,(void**)(&ozz)); 304 PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->constants[1]->ij,(void**)(&zoz)); 305 if (jac->constants[2]) { 306 PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->constants[2]->ij,(void**)(&zzo)); 307 } 308 PetscStackCallStandard(HYPRE_AMSSetEdgeConstantVectors,jac->hsolver,ozz,zoz,zzo); 309 } 310 if (jac->coords[0]) { 311 HYPRE_ParVector coords[3]; 312 coords[0] = NULL; 313 coords[1] = NULL; 314 coords[2] = NULL; 315 if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->coords[0]->ij,(void**)(&coords[0])); 316 if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->coords[1]->ij,(void**)(&coords[1])); 317 if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->coords[2]->ij,(void**)(&coords[2])); 318 PetscStackCallStandard(HYPRE_AMSSetCoordinateVectors,jac->hsolver,coords[0],coords[1],coords[2]); 319 } 320 PetscCheck(jac->G,PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE AMS preconditioner needs the discrete gradient operator via PCHYPRESetDiscreteGradient"); 321 hm = (Mat_HYPRE*)(jac->G->data); 322 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&parcsr)); 323 PetscStackCallStandard(HYPRE_AMSSetDiscreteGradient,jac->hsolver,parcsr); 324 if (jac->alpha_Poisson) { 325 hm = (Mat_HYPRE*)(jac->alpha_Poisson->data); 326 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&parcsr)); 327 PetscStackCallStandard(HYPRE_AMSSetAlphaPoissonMatrix,jac->hsolver,parcsr); 328 } 329 if (jac->ams_beta_is_zero) { 330 PetscStackCallStandard(HYPRE_AMSSetBetaPoissonMatrix,jac->hsolver,NULL); 331 } else if (jac->beta_Poisson) { 332 hm = (Mat_HYPRE*)(jac->beta_Poisson->data); 333 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&parcsr)); 334 PetscStackCallStandard(HYPRE_AMSSetBetaPoissonMatrix,jac->hsolver,parcsr); 335 } 336 if (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1])) { 337 PetscInt i; 338 HYPRE_ParCSRMatrix nd_parcsrfull, nd_parcsr[3]; 339 if (jac->ND_PiFull) { 340 hm = (Mat_HYPRE*)(jac->ND_PiFull->data); 341 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&nd_parcsrfull)); 342 } else { 343 nd_parcsrfull = NULL; 344 } 345 for (i=0;i<3;++i) { 346 if (jac->ND_Pi[i]) { 347 hm = (Mat_HYPRE*)(jac->ND_Pi[i]->data); 348 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&nd_parcsr[i])); 349 } else { 350 nd_parcsr[i] = NULL; 351 } 352 } 353 PetscStackCallStandard(HYPRE_AMSSetInterpolations,jac->hsolver,nd_parcsrfull,nd_parcsr[0],nd_parcsr[1],nd_parcsr[2]); 354 } 355 } 356 /* special case for ADS */ 357 if (jac->setup == HYPRE_ADSSetup) { 358 Mat_HYPRE *hm; 359 HYPRE_ParCSRMatrix parcsr; 360 if (!jac->coords[0] && !((jac->RT_PiFull || (jac->RT_Pi[0] && jac->RT_Pi[1])) && (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1])))) { 361 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs either the coordinate vectors via PCSetCoordinates() or the interpolation matrices via PCHYPRESetInterpolations"); 362 } 363 else PetscCheckFalse(!jac->coords[1] || !jac->coords[2],PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner has been designed for three dimensional problems! For two dimensional problems, use HYPRE AMS instead"); 364 PetscCheck(jac->G,PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs the discrete gradient operator via PCHYPRESetDiscreteGradient"); 365 PetscCheck(jac->C,PetscObjectComm((PetscObject)pc),PETSC_ERR_USER,"HYPRE ADS preconditioner needs the discrete curl operator via PCHYPRESetDiscreteGradient"); 366 if (jac->coords[0]) { 367 HYPRE_ParVector coords[3]; 368 coords[0] = NULL; 369 coords[1] = NULL; 370 coords[2] = NULL; 371 if (jac->coords[0]) PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->coords[0]->ij,(void**)(&coords[0])); 372 if (jac->coords[1]) PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->coords[1]->ij,(void**)(&coords[1])); 373 if (jac->coords[2]) PetscStackCallStandard(HYPRE_IJVectorGetObject,jac->coords[2]->ij,(void**)(&coords[2])); 374 PetscStackCallStandard(HYPRE_ADSSetCoordinateVectors,jac->hsolver,coords[0],coords[1],coords[2]); 375 } 376 hm = (Mat_HYPRE*)(jac->G->data); 377 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&parcsr)); 378 PetscStackCallStandard(HYPRE_ADSSetDiscreteGradient,jac->hsolver,parcsr); 379 hm = (Mat_HYPRE*)(jac->C->data); 380 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&parcsr)); 381 PetscStackCallStandard(HYPRE_ADSSetDiscreteCurl,jac->hsolver,parcsr); 382 if ((jac->RT_PiFull || (jac->RT_Pi[0] && jac->RT_Pi[1])) && (jac->ND_PiFull || (jac->ND_Pi[0] && jac->ND_Pi[1]))) { 383 PetscInt i; 384 HYPRE_ParCSRMatrix rt_parcsrfull, rt_parcsr[3]; 385 HYPRE_ParCSRMatrix nd_parcsrfull, nd_parcsr[3]; 386 if (jac->RT_PiFull) { 387 hm = (Mat_HYPRE*)(jac->RT_PiFull->data); 388 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&rt_parcsrfull)); 389 } else { 390 rt_parcsrfull = NULL; 391 } 392 for (i=0;i<3;++i) { 393 if (jac->RT_Pi[i]) { 394 hm = (Mat_HYPRE*)(jac->RT_Pi[i]->data); 395 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&rt_parcsr[i])); 396 } else { 397 rt_parcsr[i] = NULL; 398 } 399 } 400 if (jac->ND_PiFull) { 401 hm = (Mat_HYPRE*)(jac->ND_PiFull->data); 402 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&nd_parcsrfull)); 403 } else { 404 nd_parcsrfull = NULL; 405 } 406 for (i=0;i<3;++i) { 407 if (jac->ND_Pi[i]) { 408 hm = (Mat_HYPRE*)(jac->ND_Pi[i]->data); 409 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hm->ij,(void**)(&nd_parcsr[i])); 410 } else { 411 nd_parcsr[i] = NULL; 412 } 413 } 414 PetscStackCallStandard(HYPRE_ADSSetInterpolations,jac->hsolver,rt_parcsrfull,rt_parcsr[0],rt_parcsr[1],rt_parcsr[2],nd_parcsrfull,nd_parcsr[0],nd_parcsr[1],nd_parcsr[2]); 415 } 416 } 417 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hjac->ij,(void**)&hmat); 418 PetscStackCallStandard(HYPRE_IJVectorGetObject,hjac->b->ij,(void**)&bv); 419 PetscStackCallStandard(HYPRE_IJVectorGetObject,hjac->x->ij,(void**)&xv); 420 PetscStackCallStandard(jac->setup,jac->hsolver,hmat,bv,xv); 421 PetscFunctionReturn(0); 422 } 423 424 static PetscErrorCode PCApply_HYPRE(PC pc,Vec b,Vec x) 425 { 426 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 427 Mat_HYPRE *hjac = (Mat_HYPRE*)(jac->hpmat->data); 428 HYPRE_ParCSRMatrix hmat; 429 HYPRE_ParVector jbv,jxv; 430 431 PetscFunctionBegin; 432 CHKERRQ(PetscCitationsRegister(hypreCitation,&cite)); 433 if (!jac->applyrichardson) CHKERRQ(VecSet(x,0.0)); 434 CHKERRQ(VecHYPRE_IJVectorPushVecRead(hjac->b,b)); 435 if (jac->applyrichardson) CHKERRQ(VecHYPRE_IJVectorPushVec(hjac->x,x)); 436 else CHKERRQ(VecHYPRE_IJVectorPushVecWrite(hjac->x,x)); 437 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hjac->ij,(void**)&hmat); 438 PetscStackCallStandard(HYPRE_IJVectorGetObject,hjac->b->ij,(void**)&jbv); 439 PetscStackCallStandard(HYPRE_IJVectorGetObject,hjac->x->ij,(void**)&jxv); 440 PetscStackCall("Hypre solve",do { 441 HYPRE_Int hierr = (*jac->solve)(jac->hsolver,hmat,jbv,jxv); 442 if (hierr) { 443 PetscCheck(hierr == HYPRE_ERROR_CONV,PETSC_COMM_SELF,PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",(int)hierr); 444 hypre__global_error = 0; 445 } 446 } while (0)); 447 448 if (jac->setup == HYPRE_AMSSetup && jac->ams_beta_is_zero_part) { 449 PetscStackCallStandard(HYPRE_AMSProjectOutGradients,jac->hsolver,jxv); 450 } 451 CHKERRQ(VecHYPRE_IJVectorPopVec(hjac->x)); 452 CHKERRQ(VecHYPRE_IJVectorPopVec(hjac->b)); 453 PetscFunctionReturn(0); 454 } 455 456 static PetscErrorCode PCReset_HYPRE(PC pc) 457 { 458 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 459 460 PetscFunctionBegin; 461 CHKERRQ(MatDestroy(&jac->hpmat)); 462 CHKERRQ(MatDestroy(&jac->G)); 463 CHKERRQ(MatDestroy(&jac->C)); 464 CHKERRQ(MatDestroy(&jac->alpha_Poisson)); 465 CHKERRQ(MatDestroy(&jac->beta_Poisson)); 466 CHKERRQ(MatDestroy(&jac->RT_PiFull)); 467 CHKERRQ(MatDestroy(&jac->RT_Pi[0])); 468 CHKERRQ(MatDestroy(&jac->RT_Pi[1])); 469 CHKERRQ(MatDestroy(&jac->RT_Pi[2])); 470 CHKERRQ(MatDestroy(&jac->ND_PiFull)); 471 CHKERRQ(MatDestroy(&jac->ND_Pi[0])); 472 CHKERRQ(MatDestroy(&jac->ND_Pi[1])); 473 CHKERRQ(MatDestroy(&jac->ND_Pi[2])); 474 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->coords[0])); 475 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->coords[1])); 476 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->coords[2])); 477 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->constants[0])); 478 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->constants[1])); 479 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->constants[2])); 480 CHKERRQ(PCHYPREResetNearNullSpace_Private(pc)); 481 jac->ams_beta_is_zero = PETSC_FALSE; 482 jac->dim = 0; 483 PetscFunctionReturn(0); 484 } 485 486 static PetscErrorCode PCDestroy_HYPRE(PC pc) 487 { 488 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 489 490 PetscFunctionBegin; 491 CHKERRQ(PCReset_HYPRE(pc)); 492 if (jac->destroy) PetscStackCallStandard(jac->destroy,jac->hsolver); 493 CHKERRQ(PetscFree(jac->hypre_type)); 494 #if PETSC_PKG_HYPRE_VERSION_GE(2,23,0) 495 CHKERRQ(PetscFree(jac->spgemm_type)); 496 #endif 497 if (jac->comm_hypre != MPI_COMM_NULL) CHKERRQ(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc),&jac->comm_hypre)); 498 CHKERRQ(PetscFree(pc->data)); 499 500 CHKERRQ(PetscObjectChangeTypeName((PetscObject)pc,0)); 501 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetType_C",NULL)); 502 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPREGetType_C",NULL)); 503 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetCoordinates_C",NULL)); 504 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteGradient_C",NULL)); 505 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteCurl_C",NULL)); 506 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetInterpolations_C",NULL)); 507 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetConstantEdgeVectors_C",NULL)); 508 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetPoissonMatrix_C",NULL)); 509 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCGetInterpolations_C",NULL)); 510 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCGetCoarseOperators_C",NULL)); 511 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCMGGalerkinSetMatProductAlgorithm_C",NULL)); 512 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCMGGalerkinGetMatProductAlgorithm_C",NULL)); 513 514 PetscFunctionReturn(0); 515 } 516 517 /* --------------------------------------------------------------------------------------------*/ 518 static PetscErrorCode PCSetFromOptions_HYPRE_Pilut(PetscOptionItems *PetscOptionsObject,PC pc) 519 { 520 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 521 PetscBool flag; 522 523 PetscFunctionBegin; 524 CHKERRQ(PetscOptionsHead(PetscOptionsObject,"HYPRE Pilut Options")); 525 CHKERRQ(PetscOptionsInt("-pc_hypre_pilut_maxiter","Number of iterations","None",jac->maxiter,&jac->maxiter,&flag)); 526 if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetMaxIter,jac->hsolver,jac->maxiter); 527 CHKERRQ(PetscOptionsReal("-pc_hypre_pilut_tol","Drop tolerance","None",jac->tol,&jac->tol,&flag)); 528 if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetDropTolerance,jac->hsolver,jac->tol); 529 CHKERRQ(PetscOptionsInt("-pc_hypre_pilut_factorrowsize","FactorRowSize","None",jac->factorrowsize,&jac->factorrowsize,&flag)); 530 if (flag) PetscStackCallStandard(HYPRE_ParCSRPilutSetFactorRowSize,jac->hsolver,jac->factorrowsize); 531 CHKERRQ(PetscOptionsTail()); 532 PetscFunctionReturn(0); 533 } 534 535 static PetscErrorCode PCView_HYPRE_Pilut(PC pc,PetscViewer viewer) 536 { 537 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 538 PetscBool iascii; 539 540 PetscFunctionBegin; 541 CHKERRQ(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii)); 542 if (iascii) { 543 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE Pilut preconditioning\n")); 544 if (jac->maxiter != PETSC_DEFAULT) { 545 CHKERRQ(PetscViewerASCIIPrintf(viewer," maximum number of iterations %d\n",jac->maxiter)); 546 } else { 547 CHKERRQ(PetscViewerASCIIPrintf(viewer," default maximum number of iterations \n")); 548 } 549 if (jac->tol != PETSC_DEFAULT) { 550 CHKERRQ(PetscViewerASCIIPrintf(viewer," drop tolerance %g\n",(double)jac->tol)); 551 } else { 552 CHKERRQ(PetscViewerASCIIPrintf(viewer," default drop tolerance \n")); 553 } 554 if (jac->factorrowsize != PETSC_DEFAULT) { 555 CHKERRQ(PetscViewerASCIIPrintf(viewer," factor row size %d\n",jac->factorrowsize)); 556 } else { 557 CHKERRQ(PetscViewerASCIIPrintf(viewer," default factor row size \n")); 558 } 559 } 560 PetscFunctionReturn(0); 561 } 562 563 /* --------------------------------------------------------------------------------------------*/ 564 static PetscErrorCode PCSetFromOptions_HYPRE_Euclid(PetscOptionItems *PetscOptionsObject,PC pc) 565 { 566 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 567 PetscBool flag,eu_bj = jac->eu_bj ? PETSC_TRUE : PETSC_FALSE; 568 569 PetscFunctionBegin; 570 CHKERRQ(PetscOptionsHead(PetscOptionsObject,"HYPRE Euclid Options")); 571 CHKERRQ(PetscOptionsInt("-pc_hypre_euclid_level","Factorization levels","None",jac->eu_level,&jac->eu_level,&flag)); 572 if (flag) PetscStackCallStandard(HYPRE_EuclidSetLevel,jac->hsolver,jac->eu_level); 573 574 CHKERRQ(PetscOptionsReal("-pc_hypre_euclid_droptolerance","Drop tolerance for ILU(k) in Euclid","None",jac->eu_droptolerance,&jac->eu_droptolerance,&flag)); 575 if (flag) { 576 PetscMPIInt size; 577 578 CHKERRMPI(MPI_Comm_size(PetscObjectComm((PetscObject)pc),&size)); 579 PetscCheckFalse(size > 1,PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"hypre's Euclid does not support a parallel drop tolerance"); 580 PetscStackCallStandard(HYPRE_EuclidSetILUT,jac->hsolver,jac->eu_droptolerance); 581 } 582 583 CHKERRQ(PetscOptionsBool("-pc_hypre_euclid_bj", "Use Block Jacobi for ILU in Euclid", "None", eu_bj,&eu_bj,&flag)); 584 if (flag) { 585 jac->eu_bj = eu_bj ? 1 : 0; 586 PetscStackCallStandard(HYPRE_EuclidSetBJ,jac->hsolver,jac->eu_bj); 587 } 588 CHKERRQ(PetscOptionsTail()); 589 PetscFunctionReturn(0); 590 } 591 592 static PetscErrorCode PCView_HYPRE_Euclid(PC pc,PetscViewer viewer) 593 { 594 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 595 PetscBool iascii; 596 597 PetscFunctionBegin; 598 CHKERRQ(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii)); 599 if (iascii) { 600 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE Euclid preconditioning\n")); 601 if (jac->eu_level != PETSC_DEFAULT) { 602 CHKERRQ(PetscViewerASCIIPrintf(viewer," factorization levels %d\n",jac->eu_level)); 603 } else { 604 CHKERRQ(PetscViewerASCIIPrintf(viewer," default factorization levels \n")); 605 } 606 CHKERRQ(PetscViewerASCIIPrintf(viewer," drop tolerance %g\n",(double)jac->eu_droptolerance)); 607 CHKERRQ(PetscViewerASCIIPrintf(viewer," use Block-Jacobi? %D\n",jac->eu_bj)); 608 } 609 PetscFunctionReturn(0); 610 } 611 612 /* --------------------------------------------------------------------------------------------*/ 613 614 static PetscErrorCode PCApplyTranspose_HYPRE_BoomerAMG(PC pc,Vec b,Vec x) 615 { 616 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 617 Mat_HYPRE *hjac = (Mat_HYPRE*)(jac->hpmat->data); 618 HYPRE_ParCSRMatrix hmat; 619 HYPRE_ParVector jbv,jxv; 620 621 PetscFunctionBegin; 622 CHKERRQ(PetscCitationsRegister(hypreCitation,&cite)); 623 CHKERRQ(VecSet(x,0.0)); 624 CHKERRQ(VecHYPRE_IJVectorPushVecRead(hjac->x,b)); 625 CHKERRQ(VecHYPRE_IJVectorPushVecWrite(hjac->b,x)); 626 627 PetscStackCallStandard(HYPRE_IJMatrixGetObject,hjac->ij,(void**)&hmat); 628 PetscStackCallStandard(HYPRE_IJVectorGetObject,hjac->b->ij,(void**)&jbv); 629 PetscStackCallStandard(HYPRE_IJVectorGetObject,hjac->x->ij,(void**)&jxv); 630 631 PetscStackCall("Hypre Transpose solve",do { 632 HYPRE_Int hierr = HYPRE_BoomerAMGSolveT(jac->hsolver,hmat,jbv,jxv); 633 if (hierr) { 634 /* error code of 1 in BoomerAMG merely means convergence not achieved */ 635 PetscCheck(hierr == 1,PETSC_COMM_SELF,PETSC_ERR_LIB,"Error in HYPRE solver, error code %d",(int)hierr); 636 hypre__global_error = 0; 637 } 638 } while (0)); 639 640 CHKERRQ(VecHYPRE_IJVectorPopVec(hjac->x)); 641 CHKERRQ(VecHYPRE_IJVectorPopVec(hjac->b)); 642 PetscFunctionReturn(0); 643 } 644 645 /* static array length */ 646 #define ALEN(a) (sizeof(a)/sizeof((a)[0])) 647 648 static PetscErrorCode PCMGGalerkinSetMatProductAlgorithm_HYPRE_BoomerAMG(PC pc,const char name[]) 649 { 650 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 651 PetscBool flag; 652 653 #if PETSC_PKG_HYPRE_VERSION_GE(2,23,0) 654 PetscFunctionBegin; 655 if (jac->spgemm_type) { 656 CHKERRQ(PetscStrcmp(jac->spgemm_type,name,&flag)); 657 PetscCheck(flag,PetscObjectComm((PetscObject)pc),PETSC_ERR_ORDER,"Cannot reset the HYPRE SpGEMM (really we can)"); 658 PetscFunctionReturn(0); 659 } else { 660 CHKERRQ(PetscStrallocpy(name, &jac->spgemm_type)); 661 } 662 CHKERRQ(PetscStrcmp("cusparse",jac->spgemm_type,&flag)); 663 if (flag) { 664 PetscStackCallStandard(HYPRE_SetSpGemmUseCusparse,1); 665 PetscFunctionReturn(0); 666 } 667 CHKERRQ(PetscStrcmp("hypre",jac->spgemm_type,&flag)); 668 if (flag) { 669 PetscStackCallStandard(HYPRE_SetSpGemmUseCusparse,0); 670 PetscFunctionReturn(0); 671 } 672 jac->spgemm_type = NULL; 673 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown HYPRE SpGEM type %s; Choices are cusparse, hypre",name); 674 #endif 675 } 676 677 static PetscErrorCode PCMGGalerkinGetMatProductAlgorithm_HYPRE_BoomerAMG(PC pc, const char *spgemm[]) 678 { 679 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 680 681 PetscFunctionBegin; 682 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 683 #if PETSC_PKG_HYPRE_VERSION_GE(2,23,0) 684 *spgemm = jac->spgemm_type; 685 #endif 686 PetscFunctionReturn(0); 687 } 688 689 static const char *HYPREBoomerAMGCycleType[] = {"","V","W"}; 690 static const char *HYPREBoomerAMGCoarsenType[] = {"CLJP","Ruge-Stueben","","modifiedRuge-Stueben","","","Falgout", "", "PMIS", "", "HMIS"}; 691 static const char *HYPREBoomerAMGMeasureType[] = {"local","global"}; 692 /* The following corresponds to HYPRE_BoomerAMGSetRelaxType which has many missing numbers in the enum */ 693 static const char *HYPREBoomerAMGSmoothType[] = {"Schwarz-smoothers","Pilut","ParaSails","Euclid"}; 694 static const char *HYPREBoomerAMGRelaxType[] = {"Jacobi","sequential-Gauss-Seidel","seqboundary-Gauss-Seidel","SOR/Jacobi","backward-SOR/Jacobi", 695 "" /* [5] hybrid chaotic Gauss-Seidel (works only with OpenMP) */,"symmetric-SOR/Jacobi", 696 "" /* 7 */,"l1scaled-SOR/Jacobi","Gaussian-elimination", 697 "" /* 10 */, "" /* 11 */, "" /* 12 */, "l1-Gauss-Seidel" /* nonsymmetric */, "backward-l1-Gauss-Seidel" /* nonsymmetric */, 698 "CG" /* non-stationary */,"Chebyshev","FCF-Jacobi","l1scaled-Jacobi"}; 699 static const char *HYPREBoomerAMGInterpType[] = {"classical", "", "", "direct", "multipass", "multipass-wts", "ext+i", 700 "ext+i-cc", "standard", "standard-wts", "block", "block-wtd", "FF", "FF1", 701 "ext", "ad-wts", "ext-mm", "ext+i-mm", "ext+e-mm"}; 702 static PetscErrorCode PCSetFromOptions_HYPRE_BoomerAMG(PetscOptionItems *PetscOptionsObject,PC pc) 703 { 704 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 705 PetscInt bs,n,indx,level; 706 PetscBool flg, tmp_truth; 707 double tmpdbl, twodbl[2]; 708 const char *symtlist[] = {"nonsymmetric","SPD","nonsymmetric,SPD"}; 709 const char *PCHYPRESpgemmTypes[] = {"cusparse","hypre"}; 710 711 PetscFunctionBegin; 712 CHKERRQ(PetscOptionsHead(PetscOptionsObject,"HYPRE BoomerAMG Options")); 713 CHKERRQ(PetscOptionsEList("-pc_hypre_boomeramg_cycle_type","Cycle type","None",HYPREBoomerAMGCycleType+1,2,HYPREBoomerAMGCycleType[jac->cycletype],&indx,&flg)); 714 if (flg) { 715 jac->cycletype = indx+1; 716 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleType,jac->hsolver,jac->cycletype); 717 } 718 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_max_levels","Number of levels (of grids) allowed","None",jac->maxlevels,&jac->maxlevels,&flg)); 719 if (flg) { 720 PetscCheckFalse(jac->maxlevels < 2,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of levels %d must be at least two",jac->maxlevels); 721 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxLevels,jac->hsolver,jac->maxlevels); 722 } 723 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_max_iter","Maximum iterations used PER hypre call","None",jac->maxiter,&jac->maxiter,&flg)); 724 if (flg) { 725 PetscCheckFalse(jac->maxiter < 1,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of iterations %d must be at least one",jac->maxiter); 726 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,jac->hsolver,jac->maxiter); 727 } 728 CHKERRQ(PetscOptionsReal("-pc_hypre_boomeramg_tol","Convergence tolerance PER hypre call (0.0 = use a fixed number of iterations)","None",jac->tol,&jac->tol,&flg)); 729 if (flg) { 730 PetscCheckFalse(jac->tol < 0.0,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Tolerance %g must be greater than or equal to zero",(double)jac->tol); 731 PetscStackCallStandard(HYPRE_BoomerAMGSetTol,jac->hsolver,jac->tol); 732 } 733 bs = 1; 734 if (pc->pmat) { 735 CHKERRQ(MatGetBlockSize(pc->pmat,&bs)); 736 } 737 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_numfunctions","Number of functions","HYPRE_BoomerAMGSetNumFunctions",bs,&bs,&flg)); 738 if (flg) { 739 PetscStackCallStandard(HYPRE_BoomerAMGSetNumFunctions,jac->hsolver,bs); 740 } 741 742 CHKERRQ(PetscOptionsReal("-pc_hypre_boomeramg_truncfactor","Truncation factor for interpolation (0=no truncation)","None",jac->truncfactor,&jac->truncfactor,&flg)); 743 if (flg) { 744 PetscCheckFalse(jac->truncfactor < 0.0,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Truncation factor %g must be great than or equal zero",(double)jac->truncfactor); 745 PetscStackCallStandard(HYPRE_BoomerAMGSetTruncFactor,jac->hsolver,jac->truncfactor); 746 } 747 748 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_P_max","Max elements per row for interpolation operator (0=unlimited)","None",jac->pmax,&jac->pmax,&flg)); 749 if (flg) { 750 PetscCheckFalse(jac->pmax < 0,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"P_max %D must be greater than or equal to zero",jac->pmax); 751 PetscStackCallStandard(HYPRE_BoomerAMGSetPMaxElmts,jac->hsolver,jac->pmax); 752 } 753 754 CHKERRQ(PetscOptionsRangeInt("-pc_hypre_boomeramg_agg_nl","Number of levels of aggressive coarsening","None",jac->agg_nl,&jac->agg_nl,&flg,0,jac->maxlevels)); 755 if (flg) PetscStackCallStandard(HYPRE_BoomerAMGSetAggNumLevels,jac->hsolver,jac->agg_nl); 756 757 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_agg_num_paths","Number of paths for aggressive coarsening","None",jac->agg_num_paths,&jac->agg_num_paths,&flg)); 758 if (flg) { 759 PetscCheckFalse(jac->agg_num_paths < 1,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Number of paths %D must be greater than or equal to 1",jac->agg_num_paths); 760 PetscStackCallStandard(HYPRE_BoomerAMGSetNumPaths,jac->hsolver,jac->agg_num_paths); 761 } 762 763 CHKERRQ(PetscOptionsReal("-pc_hypre_boomeramg_strong_threshold","Threshold for being strongly connected","None",jac->strongthreshold,&jac->strongthreshold,&flg)); 764 if (flg) { 765 PetscCheckFalse(jac->strongthreshold < 0.0,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Strong threshold %g must be great than or equal zero",(double)jac->strongthreshold); 766 PetscStackCallStandard(HYPRE_BoomerAMGSetStrongThreshold,jac->hsolver,jac->strongthreshold); 767 } 768 CHKERRQ(PetscOptionsReal("-pc_hypre_boomeramg_max_row_sum","Maximum row sum","None",jac->maxrowsum,&jac->maxrowsum,&flg)); 769 if (flg) { 770 PetscCheckFalse(jac->maxrowsum < 0.0,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Maximum row sum %g must be greater than zero",(double)jac->maxrowsum); 771 PetscCheckFalse(jac->maxrowsum > 1.0,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Maximum row sum %g must be less than or equal one",(double)jac->maxrowsum); 772 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxRowSum,jac->hsolver,jac->maxrowsum); 773 } 774 775 /* Grid sweeps */ 776 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_all","Number of sweeps for the up and down grid levels","None",jac->gridsweeps[0],&indx,&flg)); 777 if (flg) { 778 PetscStackCallStandard(HYPRE_BoomerAMGSetNumSweeps,jac->hsolver,indx); 779 /* modify the jac structure so we can view the updated options with PC_View */ 780 jac->gridsweeps[0] = indx; 781 jac->gridsweeps[1] = indx; 782 /*defaults coarse to 1 */ 783 jac->gridsweeps[2] = 1; 784 } 785 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_nodal_coarsen","Use a nodal based coarsening 1-6","HYPRE_BoomerAMGSetNodal",jac->nodal_coarsening,&jac->nodal_coarsening,&flg)); 786 if (flg) { 787 PetscStackCallStandard(HYPRE_BoomerAMGSetNodal,jac->hsolver,jac->nodal_coarsening); 788 } 789 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_nodal_coarsen_diag","Diagonal in strength matrix for nodal based coarsening 0-2","HYPRE_BoomerAMGSetNodalDiag",jac->nodal_coarsening_diag,&jac->nodal_coarsening_diag,&flg)); 790 if (flg) { 791 PetscStackCallStandard(HYPRE_BoomerAMGSetNodalDiag,jac->hsolver,jac->nodal_coarsening_diag); 792 } 793 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_vec_interp_variant","Variant of algorithm 1-3","HYPRE_BoomerAMGSetInterpVecVariant",jac->vec_interp_variant, &jac->vec_interp_variant,&flg)); 794 if (flg) { 795 PetscStackCallStandard(HYPRE_BoomerAMGSetInterpVecVariant,jac->hsolver,jac->vec_interp_variant); 796 } 797 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_vec_interp_qmax","Max elements per row for each Q","HYPRE_BoomerAMGSetInterpVecQMax",jac->vec_interp_qmax, &jac->vec_interp_qmax,&flg)); 798 if (flg) { 799 PetscStackCallStandard(HYPRE_BoomerAMGSetInterpVecQMax,jac->hsolver,jac->vec_interp_qmax); 800 } 801 CHKERRQ(PetscOptionsBool("-pc_hypre_boomeramg_vec_interp_smooth","Whether to smooth the interpolation vectors","HYPRE_BoomerAMGSetSmoothInterpVectors",jac->vec_interp_smooth, &jac->vec_interp_smooth,&flg)); 802 if (flg) { 803 PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothInterpVectors,jac->hsolver,jac->vec_interp_smooth); 804 } 805 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_interp_refine","Preprocess the interpolation matrix through iterative weight refinement","HYPRE_BoomerAMGSetInterpRefine",jac->interp_refine, &jac->interp_refine,&flg)); 806 if (flg) { 807 PetscStackCallStandard(HYPRE_BoomerAMGSetInterpRefine,jac->hsolver,jac->interp_refine); 808 } 809 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_down","Number of sweeps for the down cycles","None",jac->gridsweeps[0], &indx,&flg)); 810 if (flg) { 811 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,jac->hsolver,indx, 1); 812 jac->gridsweeps[0] = indx; 813 } 814 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_up","Number of sweeps for the up cycles","None",jac->gridsweeps[1],&indx,&flg)); 815 if (flg) { 816 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,jac->hsolver,indx, 2); 817 jac->gridsweeps[1] = indx; 818 } 819 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_grid_sweeps_coarse","Number of sweeps for the coarse level","None",jac->gridsweeps[2],&indx,&flg)); 820 if (flg) { 821 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleNumSweeps,jac->hsolver,indx, 3); 822 jac->gridsweeps[2] = indx; 823 } 824 825 /* Smooth type */ 826 CHKERRQ(PetscOptionsEList("-pc_hypre_boomeramg_smooth_type","Enable more complex smoothers","None",HYPREBoomerAMGSmoothType,ALEN(HYPREBoomerAMGSmoothType),HYPREBoomerAMGSmoothType[0],&indx,&flg)); 827 if (flg) { 828 jac->smoothtype = indx; 829 PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothType,jac->hsolver,indx+6); 830 jac->smoothnumlevels = 25; 831 PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothNumLevels,jac->hsolver,25); 832 } 833 834 /* Number of smoothing levels */ 835 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_smooth_num_levels","Number of levels on which more complex smoothers are used","None",25,&indx,&flg)); 836 if (flg && (jac->smoothtype != -1)) { 837 jac->smoothnumlevels = indx; 838 PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothNumLevels,jac->hsolver,indx); 839 } 840 841 /* Number of levels for ILU(k) for Euclid */ 842 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_eu_level","Number of levels for ILU(k) in Euclid smoother","None",0,&indx,&flg)); 843 if (flg && (jac->smoothtype == 3)) { 844 jac->eu_level = indx; 845 PetscStackCallStandard(HYPRE_BoomerAMGSetEuLevel,jac->hsolver,indx); 846 } 847 848 /* Filter for ILU(k) for Euclid */ 849 double droptolerance; 850 CHKERRQ(PetscOptionsReal("-pc_hypre_boomeramg_eu_droptolerance","Drop tolerance for ILU(k) in Euclid smoother","None",0,&droptolerance,&flg)); 851 if (flg && (jac->smoothtype == 3)) { 852 jac->eu_droptolerance = droptolerance; 853 PetscStackCallStandard(HYPRE_BoomerAMGSetEuLevel,jac->hsolver,droptolerance); 854 } 855 856 /* Use Block Jacobi ILUT for Euclid */ 857 CHKERRQ(PetscOptionsBool("-pc_hypre_boomeramg_eu_bj", "Use Block Jacobi for ILU in Euclid smoother?", "None", PETSC_FALSE, &tmp_truth, &flg)); 858 if (flg && (jac->smoothtype == 3)) { 859 jac->eu_bj = tmp_truth; 860 PetscStackCallStandard(HYPRE_BoomerAMGSetEuBJ,jac->hsolver,jac->eu_bj); 861 } 862 863 /* Relax type */ 864 CHKERRQ(PetscOptionsEList("-pc_hypre_boomeramg_relax_type_all","Relax type for the up and down cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg)); 865 if (flg) { 866 jac->relaxtype[0] = jac->relaxtype[1] = indx; 867 PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxType,jac->hsolver, indx); 868 /* by default, coarse type set to 9 */ 869 jac->relaxtype[2] = 9; 870 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,jac->hsolver, 9, 3); 871 } 872 CHKERRQ(PetscOptionsEList("-pc_hypre_boomeramg_relax_type_down","Relax type for the down cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg)); 873 if (flg) { 874 jac->relaxtype[0] = indx; 875 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,jac->hsolver, indx, 1); 876 } 877 CHKERRQ(PetscOptionsEList("-pc_hypre_boomeramg_relax_type_up","Relax type for the up cycles","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[6],&indx,&flg)); 878 if (flg) { 879 jac->relaxtype[1] = indx; 880 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,jac->hsolver, indx, 2); 881 } 882 CHKERRQ(PetscOptionsEList("-pc_hypre_boomeramg_relax_type_coarse","Relax type on coarse grid","None",HYPREBoomerAMGRelaxType,ALEN(HYPREBoomerAMGRelaxType),HYPREBoomerAMGRelaxType[9],&indx,&flg)); 883 if (flg) { 884 jac->relaxtype[2] = indx; 885 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleRelaxType,jac->hsolver, indx, 3); 886 } 887 888 /* Relaxation Weight */ 889 CHKERRQ(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)); 890 if (flg) { 891 PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxWt,jac->hsolver,tmpdbl); 892 jac->relaxweight = tmpdbl; 893 } 894 895 n = 2; 896 twodbl[0] = twodbl[1] = 1.0; 897 CHKERRQ(PetscOptionsRealArray("-pc_hypre_boomeramg_relax_weight_level","Set the relaxation weight for a particular level (weight,level)","None",twodbl, &n, &flg)); 898 if (flg) { 899 if (n == 2) { 900 indx = (int)PetscAbsReal(twodbl[1]); 901 PetscStackCallStandard(HYPRE_BoomerAMGSetLevelRelaxWt,jac->hsolver,twodbl[0],indx); 902 } else SETERRQ(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); 903 } 904 905 /* Outer relaxation Weight */ 906 CHKERRQ(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)); 907 if (flg) { 908 PetscStackCallStandard(HYPRE_BoomerAMGSetOuterWt,jac->hsolver, tmpdbl); 909 jac->outerrelaxweight = tmpdbl; 910 } 911 912 n = 2; 913 twodbl[0] = twodbl[1] = 1.0; 914 CHKERRQ(PetscOptionsRealArray("-pc_hypre_boomeramg_outer_relax_weight_level","Set the outer relaxation weight for a particular level (weight,level)","None",twodbl, &n, &flg)); 915 if (flg) { 916 if (n == 2) { 917 indx = (int)PetscAbsReal(twodbl[1]); 918 PetscStackCallStandard(HYPRE_BoomerAMGSetLevelOuterWt,jac->hsolver, twodbl[0], indx); 919 } else SETERRQ(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); 920 } 921 922 /* the Relax Order */ 923 CHKERRQ(PetscOptionsBool("-pc_hypre_boomeramg_no_CF", "Do not use CF-relaxation", "None", PETSC_FALSE, &tmp_truth, &flg)); 924 925 if (flg && tmp_truth) { 926 jac->relaxorder = 0; 927 PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxOrder,jac->hsolver, jac->relaxorder); 928 } 929 CHKERRQ(PetscOptionsEList("-pc_hypre_boomeramg_measure_type","Measure type","None",HYPREBoomerAMGMeasureType,ALEN(HYPREBoomerAMGMeasureType),HYPREBoomerAMGMeasureType[0],&indx,&flg)); 930 if (flg) { 931 jac->measuretype = indx; 932 PetscStackCallStandard(HYPRE_BoomerAMGSetMeasureType,jac->hsolver,jac->measuretype); 933 } 934 /* update list length 3/07 */ 935 CHKERRQ(PetscOptionsEList("-pc_hypre_boomeramg_coarsen_type","Coarsen type","None",HYPREBoomerAMGCoarsenType,ALEN(HYPREBoomerAMGCoarsenType),HYPREBoomerAMGCoarsenType[6],&indx,&flg)); 936 if (flg) { 937 jac->coarsentype = indx; 938 PetscStackCallStandard(HYPRE_BoomerAMGSetCoarsenType,jac->hsolver,jac->coarsentype); 939 } 940 941 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_max_coarse_size", "Maximum size of coarsest grid", "None", jac->maxc, &jac->maxc, &flg)); 942 if (flg) { 943 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxCoarseSize,jac->hsolver, jac->maxc); 944 } 945 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_min_coarse_size", "Minimum size of coarsest grid", "None", jac->minc, &jac->minc, &flg)); 946 if (flg) { 947 PetscStackCallStandard(HYPRE_BoomerAMGSetMinCoarseSize,jac->hsolver, jac->minc); 948 } 949 #if PETSC_PKG_HYPRE_VERSION_GE(2,23,0) 950 // global parameter but is closely associated with BoomerAMG 951 CHKERRQ(PetscOptionsEList("-pc_mg_galerkin_mat_product_algorithm","Type of SpGEMM to use in hypre (only for now)","PCMGGalerkinSetMatProductAlgorithm",PCHYPRESpgemmTypes,ALEN(PCHYPRESpgemmTypes),PCHYPRESpgemmTypes[0],&indx,&flg)); 952 if (!flg) indx = 0; 953 CHKERRQ(PCMGGalerkinSetMatProductAlgorithm_HYPRE_BoomerAMG(pc,PCHYPRESpgemmTypes[indx])); 954 #endif 955 /* AIR */ 956 #if PETSC_PKG_HYPRE_VERSION_GE(2,18,0) 957 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_restriction_type", "Type of AIR method (distance 1 or 2, 0 means no AIR)", "None", jac->Rtype, &jac->Rtype, NULL)); 958 PetscStackCallStandard(HYPRE_BoomerAMGSetRestriction,jac->hsolver,jac->Rtype); 959 if (jac->Rtype) { 960 jac->interptype = 100; /* no way we can pass this with strings... Set it as default as in MFEM, then users can still customize it back to a different one */ 961 962 CHKERRQ(PetscOptionsReal("-pc_hypre_boomeramg_strongthresholdR","Threshold for R","None",jac->Rstrongthreshold,&jac->Rstrongthreshold,NULL)); 963 PetscStackCallStandard(HYPRE_BoomerAMGSetStrongThresholdR,jac->hsolver,jac->Rstrongthreshold); 964 965 CHKERRQ(PetscOptionsReal("-pc_hypre_boomeramg_filterthresholdR","Filter threshold for R","None",jac->Rfilterthreshold,&jac->Rfilterthreshold,NULL)); 966 PetscStackCallStandard(HYPRE_BoomerAMGSetFilterThresholdR,jac->hsolver,jac->Rfilterthreshold); 967 968 CHKERRQ(PetscOptionsReal("-pc_hypre_boomeramg_Adroptol","Defines the drop tolerance for the A-matrices from the 2nd level of AMG","None",jac->Adroptol,&jac->Adroptol,NULL)); 969 PetscStackCallStandard(HYPRE_BoomerAMGSetADropTol,jac->hsolver,jac->Adroptol); 970 971 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_Adroptype","Drops the entries that are not on the diagonal and smaller than its row norm: type 1: 1-norm, 2: 2-norm, -1: infinity norm","None",jac->Adroptype,&jac->Adroptype,NULL)); 972 PetscStackCallStandard(HYPRE_BoomerAMGSetADropType,jac->hsolver,jac->Adroptype); 973 } 974 #endif 975 976 #if PETSC_PKG_HYPRE_VERSION_LE(9,9,9) 977 PetscCheckFalse(jac->Rtype && jac->agg_nl,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"-pc_hypre_boomeramg_restriction_type (%D) and -pc_hypre_boomeramg_agg_nl (%D)",jac->Rtype,jac->agg_nl); 978 #endif 979 980 /* new 3/07 */ 981 CHKERRQ(PetscOptionsEList("-pc_hypre_boomeramg_interp_type","Interpolation type","None",HYPREBoomerAMGInterpType,ALEN(HYPREBoomerAMGInterpType),HYPREBoomerAMGInterpType[0],&indx,&flg)); 982 if (flg || jac->Rtype) { 983 if (flg) jac->interptype = indx; 984 PetscStackCallStandard(HYPRE_BoomerAMGSetInterpType,jac->hsolver,jac->interptype); 985 } 986 987 CHKERRQ(PetscOptionsName("-pc_hypre_boomeramg_print_statistics","Print statistics","None",&flg)); 988 if (flg) { 989 level = 3; 990 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_print_statistics","Print statistics","None",level,&level,NULL)); 991 992 jac->printstatistics = PETSC_TRUE; 993 PetscStackCallStandard(HYPRE_BoomerAMGSetPrintLevel,jac->hsolver,level); 994 } 995 996 CHKERRQ(PetscOptionsName("-pc_hypre_boomeramg_print_debug","Print debug information","None",&flg)); 997 if (flg) { 998 level = 3; 999 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_print_debug","Print debug information","None",level,&level,NULL)); 1000 1001 jac->printstatistics = PETSC_TRUE; 1002 PetscStackCallStandard(HYPRE_BoomerAMGSetDebugFlag,jac->hsolver,level); 1003 } 1004 1005 CHKERRQ(PetscOptionsBool("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None", PETSC_FALSE, &tmp_truth, &flg)); 1006 if (flg && tmp_truth) { 1007 PetscInt tmp_int; 1008 CHKERRQ(PetscOptionsInt("-pc_hypre_boomeramg_nodal_relaxation", "Nodal relaxation via Schwarz", "None",jac->nodal_relax_levels,&tmp_int,&flg)); 1009 if (flg) jac->nodal_relax_levels = tmp_int; 1010 PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothType,jac->hsolver,6); 1011 PetscStackCallStandard(HYPRE_BoomerAMGSetDomainType,jac->hsolver,1); 1012 PetscStackCallStandard(HYPRE_BoomerAMGSetOverlap,jac->hsolver,0); 1013 PetscStackCallStandard(HYPRE_BoomerAMGSetSmoothNumLevels,jac->hsolver,jac->nodal_relax_levels); 1014 } 1015 1016 CHKERRQ(PetscOptionsBool("-pc_hypre_boomeramg_keeptranspose", "Avoid transpose matvecs in preconditioner application", "None", jac->keeptranspose, &jac->keeptranspose, NULL)); 1017 PetscStackCallStandard(HYPRE_BoomerAMGSetKeepTranspose,jac->hsolver,jac->keeptranspose ? 1 : 0); 1018 1019 /* options for ParaSails solvers */ 1020 CHKERRQ(PetscOptionsEList("-pc_hypre_boomeramg_parasails_sym","Symmetry of matrix and preconditioner","None",symtlist,ALEN(symtlist),symtlist[0],&indx,&flg)); 1021 if (flg) { 1022 jac->symt = indx; 1023 PetscStackCallStandard(HYPRE_BoomerAMGSetSym,jac->hsolver,jac->symt); 1024 } 1025 1026 CHKERRQ(PetscOptionsTail()); 1027 PetscFunctionReturn(0); 1028 } 1029 1030 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) 1031 { 1032 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1033 HYPRE_Int oits; 1034 1035 PetscFunctionBegin; 1036 CHKERRQ(PetscCitationsRegister(hypreCitation,&cite)); 1037 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,jac->hsolver,its*jac->maxiter); 1038 PetscStackCallStandard(HYPRE_BoomerAMGSetTol,jac->hsolver,rtol); 1039 jac->applyrichardson = PETSC_TRUE; 1040 CHKERRQ(PCApply_HYPRE(pc,b,y)); 1041 jac->applyrichardson = PETSC_FALSE; 1042 PetscStackCallStandard(HYPRE_BoomerAMGGetNumIterations,jac->hsolver,&oits); 1043 *outits = oits; 1044 if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS; 1045 else *reason = PCRICHARDSON_CONVERGED_RTOL; 1046 PetscStackCallStandard(HYPRE_BoomerAMGSetTol,jac->hsolver,jac->tol); 1047 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,jac->hsolver,jac->maxiter); 1048 PetscFunctionReturn(0); 1049 } 1050 1051 static PetscErrorCode PCView_HYPRE_BoomerAMG(PC pc,PetscViewer viewer) 1052 { 1053 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1054 PetscBool iascii; 1055 1056 PetscFunctionBegin; 1057 CHKERRQ(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii)); 1058 if (iascii) { 1059 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG preconditioning\n")); 1060 CHKERRQ(PetscViewerASCIIPrintf(viewer," Cycle type %s\n",HYPREBoomerAMGCycleType[jac->cycletype])); 1061 CHKERRQ(PetscViewerASCIIPrintf(viewer," Maximum number of levels %D\n",jac->maxlevels)); 1062 CHKERRQ(PetscViewerASCIIPrintf(viewer," Maximum number of iterations PER hypre call %D\n",jac->maxiter)); 1063 CHKERRQ(PetscViewerASCIIPrintf(viewer," Convergence tolerance PER hypre call %g\n",(double)jac->tol)); 1064 CHKERRQ(PetscViewerASCIIPrintf(viewer," Threshold for strong coupling %g\n",(double)jac->strongthreshold)); 1065 CHKERRQ(PetscViewerASCIIPrintf(viewer," Interpolation truncation factor %g\n",(double)jac->truncfactor)); 1066 CHKERRQ(PetscViewerASCIIPrintf(viewer," Interpolation: max elements per row %D\n",jac->pmax)); 1067 if (jac->interp_refine) { 1068 CHKERRQ(PetscViewerASCIIPrintf(viewer," Interpolation: number of steps of weighted refinement %D\n",jac->interp_refine)); 1069 } 1070 CHKERRQ(PetscViewerASCIIPrintf(viewer," Number of levels of aggressive coarsening %D\n",jac->agg_nl)); 1071 CHKERRQ(PetscViewerASCIIPrintf(viewer," Number of paths for aggressive coarsening %D\n",jac->agg_num_paths)); 1072 1073 CHKERRQ(PetscViewerASCIIPrintf(viewer," Maximum row sums %g\n",(double)jac->maxrowsum)); 1074 1075 CHKERRQ(PetscViewerASCIIPrintf(viewer," Sweeps down %D\n",jac->gridsweeps[0])); 1076 CHKERRQ(PetscViewerASCIIPrintf(viewer," Sweeps up %D\n",jac->gridsweeps[1])); 1077 CHKERRQ(PetscViewerASCIIPrintf(viewer," Sweeps on coarse %D\n",jac->gridsweeps[2])); 1078 1079 CHKERRQ(PetscViewerASCIIPrintf(viewer," Relax down %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[0]])); 1080 CHKERRQ(PetscViewerASCIIPrintf(viewer," Relax up %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[1]])); 1081 CHKERRQ(PetscViewerASCIIPrintf(viewer," Relax on coarse %s\n",HYPREBoomerAMGRelaxType[jac->relaxtype[2]])); 1082 1083 CHKERRQ(PetscViewerASCIIPrintf(viewer," Relax weight (all) %g\n",(double)jac->relaxweight)); 1084 CHKERRQ(PetscViewerASCIIPrintf(viewer," Outer relax weight (all) %g\n",(double)jac->outerrelaxweight)); 1085 1086 if (jac->relaxorder) { 1087 CHKERRQ(PetscViewerASCIIPrintf(viewer," Using CF-relaxation\n")); 1088 } else { 1089 CHKERRQ(PetscViewerASCIIPrintf(viewer," Not using CF-relaxation\n")); 1090 } 1091 if (jac->smoothtype!=-1) { 1092 CHKERRQ(PetscViewerASCIIPrintf(viewer," Smooth type %s\n",HYPREBoomerAMGSmoothType[jac->smoothtype])); 1093 CHKERRQ(PetscViewerASCIIPrintf(viewer," Smooth num levels %D\n",jac->smoothnumlevels)); 1094 } else { 1095 CHKERRQ(PetscViewerASCIIPrintf(viewer," Not using more complex smoothers.\n")); 1096 } 1097 if (jac->smoothtype==3) { 1098 CHKERRQ(PetscViewerASCIIPrintf(viewer," Euclid ILU(k) levels %D\n",jac->eu_level)); 1099 CHKERRQ(PetscViewerASCIIPrintf(viewer," Euclid ILU(k) drop tolerance %g\n",(double)jac->eu_droptolerance)); 1100 CHKERRQ(PetscViewerASCIIPrintf(viewer," Euclid ILU use Block-Jacobi? %D\n",jac->eu_bj)); 1101 } 1102 CHKERRQ(PetscViewerASCIIPrintf(viewer," Measure type %s\n",HYPREBoomerAMGMeasureType[jac->measuretype])); 1103 CHKERRQ(PetscViewerASCIIPrintf(viewer," Coarsen type %s\n",HYPREBoomerAMGCoarsenType[jac->coarsentype])); 1104 CHKERRQ(PetscViewerASCIIPrintf(viewer," Interpolation type %s\n",jac->interptype != 100 ? HYPREBoomerAMGInterpType[jac->interptype] : "1pt")); 1105 if (jac->nodal_coarsening) { 1106 CHKERRQ(PetscViewerASCIIPrintf(viewer," Using nodal coarsening with HYPRE_BOOMERAMGSetNodal() %D\n",jac->nodal_coarsening)); 1107 } 1108 if (jac->vec_interp_variant) { 1109 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE_BoomerAMGSetInterpVecVariant() %D\n",jac->vec_interp_variant)); 1110 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE_BoomerAMGSetInterpVecQMax() %D\n",jac->vec_interp_qmax)); 1111 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE_BoomerAMGSetSmoothInterpVectors() %d\n",jac->vec_interp_smooth)); 1112 } 1113 if (jac->nodal_relax) { 1114 CHKERRQ(PetscViewerASCIIPrintf(viewer," Using nodal relaxation via Schwarz smoothing on levels %D\n",jac->nodal_relax_levels)); 1115 } 1116 #if PETSC_PKG_HYPRE_VERSION_GE(2,23,0) 1117 CHKERRQ(PetscViewerASCIIPrintf(viewer," SpGEMM type %s\n",jac->spgemm_type)); 1118 #endif 1119 /* AIR */ 1120 if (jac->Rtype) { 1121 CHKERRQ(PetscViewerASCIIPrintf(viewer," Using approximate ideal restriction type %D\n",jac->Rtype)); 1122 CHKERRQ(PetscViewerASCIIPrintf(viewer," Threshold for R %g\n",(double)jac->Rstrongthreshold)); 1123 CHKERRQ(PetscViewerASCIIPrintf(viewer," Filter for R %g\n",(double)jac->Rfilterthreshold)); 1124 CHKERRQ(PetscViewerASCIIPrintf(viewer," A drop tolerance %g\n",(double)jac->Adroptol)); 1125 CHKERRQ(PetscViewerASCIIPrintf(viewer," A drop type %D\n",jac->Adroptype)); 1126 } 1127 } 1128 PetscFunctionReturn(0); 1129 } 1130 1131 /* --------------------------------------------------------------------------------------------*/ 1132 static PetscErrorCode PCSetFromOptions_HYPRE_ParaSails(PetscOptionItems *PetscOptionsObject,PC pc) 1133 { 1134 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1135 PetscInt indx; 1136 PetscBool flag; 1137 const char *symtlist[] = {"nonsymmetric","SPD","nonsymmetric,SPD"}; 1138 1139 PetscFunctionBegin; 1140 CHKERRQ(PetscOptionsHead(PetscOptionsObject,"HYPRE ParaSails Options")); 1141 CHKERRQ(PetscOptionsInt("-pc_hypre_parasails_nlevels","Number of number of levels","None",jac->nlevels,&jac->nlevels,0)); 1142 CHKERRQ(PetscOptionsReal("-pc_hypre_parasails_thresh","Threshold","None",jac->threshold,&jac->threshold,&flag)); 1143 if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetParams,jac->hsolver,jac->threshold,jac->nlevels); 1144 1145 CHKERRQ(PetscOptionsReal("-pc_hypre_parasails_filter","filter","None",jac->filter,&jac->filter,&flag)); 1146 if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetFilter,jac->hsolver,jac->filter); 1147 1148 CHKERRQ(PetscOptionsReal("-pc_hypre_parasails_loadbal","Load balance","None",jac->loadbal,&jac->loadbal,&flag)); 1149 if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetLoadbal,jac->hsolver,jac->loadbal); 1150 1151 CHKERRQ(PetscOptionsBool("-pc_hypre_parasails_logging","Print info to screen","None",(PetscBool)jac->logging,(PetscBool*)&jac->logging,&flag)); 1152 if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetLogging,jac->hsolver,jac->logging); 1153 1154 CHKERRQ(PetscOptionsBool("-pc_hypre_parasails_reuse","Reuse nonzero pattern in preconditioner","None",(PetscBool)jac->ruse,(PetscBool*)&jac->ruse,&flag)); 1155 if (flag) PetscStackCallStandard(HYPRE_ParaSailsSetReuse,jac->hsolver,jac->ruse); 1156 1157 CHKERRQ(PetscOptionsEList("-pc_hypre_parasails_sym","Symmetry of matrix and preconditioner","None",symtlist,ALEN(symtlist),symtlist[0],&indx,&flag)); 1158 if (flag) { 1159 jac->symt = indx; 1160 PetscStackCallStandard(HYPRE_ParaSailsSetSym,jac->hsolver,jac->symt); 1161 } 1162 1163 CHKERRQ(PetscOptionsTail()); 1164 PetscFunctionReturn(0); 1165 } 1166 1167 static PetscErrorCode PCView_HYPRE_ParaSails(PC pc,PetscViewer viewer) 1168 { 1169 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1170 PetscBool iascii; 1171 const char *symt = 0; 1172 1173 PetscFunctionBegin; 1174 CHKERRQ(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii)); 1175 if (iascii) { 1176 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE ParaSails preconditioning\n")); 1177 CHKERRQ(PetscViewerASCIIPrintf(viewer," nlevels %d\n",jac->nlevels)); 1178 CHKERRQ(PetscViewerASCIIPrintf(viewer," threshold %g\n",(double)jac->threshold)); 1179 CHKERRQ(PetscViewerASCIIPrintf(viewer," filter %g\n",(double)jac->filter)); 1180 CHKERRQ(PetscViewerASCIIPrintf(viewer," load balance %g\n",(double)jac->loadbal)); 1181 CHKERRQ(PetscViewerASCIIPrintf(viewer," reuse nonzero structure %s\n",PetscBools[jac->ruse])); 1182 CHKERRQ(PetscViewerASCIIPrintf(viewer," print info to screen %s\n",PetscBools[jac->logging])); 1183 if (!jac->symt) symt = "nonsymmetric matrix and preconditioner"; 1184 else if (jac->symt == 1) symt = "SPD matrix and preconditioner"; 1185 else if (jac->symt == 2) symt = "nonsymmetric matrix but SPD preconditioner"; 1186 else SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG,"Unknown HYPRE ParaSails symmetric option %d",jac->symt); 1187 CHKERRQ(PetscViewerASCIIPrintf(viewer," %s\n",symt)); 1188 } 1189 PetscFunctionReturn(0); 1190 } 1191 /* --------------------------------------------------------------------------------------------*/ 1192 static PetscErrorCode PCSetFromOptions_HYPRE_AMS(PetscOptionItems *PetscOptionsObject,PC pc) 1193 { 1194 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1195 PetscInt n; 1196 PetscBool flag,flag2,flag3,flag4; 1197 1198 PetscFunctionBegin; 1199 CHKERRQ(PetscOptionsHead(PetscOptionsObject,"HYPRE AMS Options")); 1200 CHKERRQ(PetscOptionsInt("-pc_hypre_ams_print_level","Debugging output level for AMS","None",jac->as_print,&jac->as_print,&flag)); 1201 if (flag) PetscStackCallStandard(HYPRE_AMSSetPrintLevel,jac->hsolver,jac->as_print); 1202 CHKERRQ(PetscOptionsInt("-pc_hypre_ams_max_iter","Maximum number of AMS multigrid iterations within PCApply","None",jac->as_max_iter,&jac->as_max_iter,&flag)); 1203 if (flag) PetscStackCallStandard(HYPRE_AMSSetMaxIter,jac->hsolver,jac->as_max_iter); 1204 CHKERRQ(PetscOptionsInt("-pc_hypre_ams_cycle_type","Cycle type for AMS multigrid","None",jac->ams_cycle_type,&jac->ams_cycle_type,&flag)); 1205 if (flag) PetscStackCallStandard(HYPRE_AMSSetCycleType,jac->hsolver,jac->ams_cycle_type); 1206 CHKERRQ(PetscOptionsReal("-pc_hypre_ams_tol","Error tolerance for AMS multigrid","None",jac->as_tol,&jac->as_tol,&flag)); 1207 if (flag) PetscStackCallStandard(HYPRE_AMSSetTol,jac->hsolver,jac->as_tol); 1208 CHKERRQ(PetscOptionsInt("-pc_hypre_ams_relax_type","Relaxation type for AMS smoother","None",jac->as_relax_type,&jac->as_relax_type,&flag)); 1209 CHKERRQ(PetscOptionsInt("-pc_hypre_ams_relax_times","Number of relaxation steps for AMS smoother","None",jac->as_relax_times,&jac->as_relax_times,&flag2)); 1210 CHKERRQ(PetscOptionsReal("-pc_hypre_ams_relax_weight","Relaxation weight for AMS smoother","None",jac->as_relax_weight,&jac->as_relax_weight,&flag3)); 1211 CHKERRQ(PetscOptionsReal("-pc_hypre_ams_omega","SSOR coefficient for AMS smoother","None",jac->as_omega,&jac->as_omega,&flag4)); 1212 if (flag || flag2 || flag3 || flag4) { 1213 PetscStackCallStandard(HYPRE_AMSSetSmoothingOptions,jac->hsolver,jac->as_relax_type, 1214 jac->as_relax_times, 1215 jac->as_relax_weight, 1216 jac->as_omega); 1217 } 1218 CHKERRQ(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)); 1219 n = 5; 1220 CHKERRQ(PetscOptionsIntArray("-pc_hypre_ams_amg_alpha_options","AMG options for vector Poisson","None",jac->as_amg_alpha_opts,&n,&flag2)); 1221 if (flag || flag2) { 1222 PetscStackCallStandard(HYPRE_AMSSetAlphaAMGOptions,jac->hsolver,jac->as_amg_alpha_opts[0], /* AMG coarsen type */ 1223 jac->as_amg_alpha_opts[1], /* AMG agg_levels */ 1224 jac->as_amg_alpha_opts[2], /* AMG relax_type */ 1225 jac->as_amg_alpha_theta, 1226 jac->as_amg_alpha_opts[3], /* AMG interp_type */ 1227 jac->as_amg_alpha_opts[4]); /* AMG Pmax */ 1228 } 1229 CHKERRQ(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)); 1230 n = 5; 1231 CHKERRQ(PetscOptionsIntArray("-pc_hypre_ams_amg_beta_options","AMG options for scalar Poisson solver","None",jac->as_amg_beta_opts,&n,&flag2)); 1232 if (flag || flag2) { 1233 PetscStackCallStandard(HYPRE_AMSSetBetaAMGOptions,jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */ 1234 jac->as_amg_beta_opts[1], /* AMG agg_levels */ 1235 jac->as_amg_beta_opts[2], /* AMG relax_type */ 1236 jac->as_amg_beta_theta, 1237 jac->as_amg_beta_opts[3], /* AMG interp_type */ 1238 jac->as_amg_beta_opts[4]); /* AMG Pmax */ 1239 } 1240 CHKERRQ(PetscOptionsInt("-pc_hypre_ams_projection_frequency","Frequency at which a projection onto the compatible subspace for problems with zero conductivity regions is performed","None",jac->ams_proj_freq,&jac->ams_proj_freq,&flag)); 1241 if (flag) { /* override HYPRE's default only if the options is used */ 1242 PetscStackCallStandard(HYPRE_AMSSetProjectionFrequency,jac->hsolver,jac->ams_proj_freq); 1243 } 1244 CHKERRQ(PetscOptionsTail()); 1245 PetscFunctionReturn(0); 1246 } 1247 1248 static PetscErrorCode PCView_HYPRE_AMS(PC pc,PetscViewer viewer) 1249 { 1250 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1251 PetscBool iascii; 1252 1253 PetscFunctionBegin; 1254 CHKERRQ(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii)); 1255 if (iascii) { 1256 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE AMS preconditioning\n")); 1257 CHKERRQ(PetscViewerASCIIPrintf(viewer," subspace iterations per application %d\n",jac->as_max_iter)); 1258 CHKERRQ(PetscViewerASCIIPrintf(viewer," subspace cycle type %d\n",jac->ams_cycle_type)); 1259 CHKERRQ(PetscViewerASCIIPrintf(viewer," subspace iteration tolerance %g\n",jac->as_tol)); 1260 CHKERRQ(PetscViewerASCIIPrintf(viewer," smoother type %d\n",jac->as_relax_type)); 1261 CHKERRQ(PetscViewerASCIIPrintf(viewer," number of smoothing steps %d\n",jac->as_relax_times)); 1262 CHKERRQ(PetscViewerASCIIPrintf(viewer," smoother weight %g\n",jac->as_relax_weight)); 1263 CHKERRQ(PetscViewerASCIIPrintf(viewer," smoother omega %g\n",jac->as_omega)); 1264 if (jac->alpha_Poisson) { 1265 CHKERRQ(PetscViewerASCIIPrintf(viewer," vector Poisson solver (passed in by user)\n")); 1266 } else { 1267 CHKERRQ(PetscViewerASCIIPrintf(viewer," vector Poisson solver (computed) \n")); 1268 } 1269 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG coarsening type %d\n",jac->as_amg_alpha_opts[0])); 1270 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG levels of aggressive coarsening %d\n",jac->as_amg_alpha_opts[1])); 1271 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG relaxation type %d\n",jac->as_amg_alpha_opts[2])); 1272 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG interpolation type %d\n",jac->as_amg_alpha_opts[3])); 1273 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG max nonzero elements in interpolation rows %d\n",jac->as_amg_alpha_opts[4])); 1274 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG strength threshold %g\n",jac->as_amg_alpha_theta)); 1275 if (!jac->ams_beta_is_zero) { 1276 if (jac->beta_Poisson) { 1277 CHKERRQ(PetscViewerASCIIPrintf(viewer," scalar Poisson solver (passed in by user)\n")); 1278 } else { 1279 CHKERRQ(PetscViewerASCIIPrintf(viewer," scalar Poisson solver (computed) \n")); 1280 } 1281 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG coarsening type %d\n",jac->as_amg_beta_opts[0])); 1282 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG levels of aggressive coarsening %d\n",jac->as_amg_beta_opts[1])); 1283 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG relaxation type %d\n",jac->as_amg_beta_opts[2])); 1284 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG interpolation type %d\n",jac->as_amg_beta_opts[3])); 1285 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG max nonzero elements in interpolation rows %d\n",jac->as_amg_beta_opts[4])); 1286 CHKERRQ(PetscViewerASCIIPrintf(viewer," boomerAMG strength threshold %g\n",jac->as_amg_beta_theta)); 1287 if (jac->ams_beta_is_zero_part) { 1288 CHKERRQ(PetscViewerASCIIPrintf(viewer," compatible subspace projection frequency %d (-1 HYPRE uses default)\n",jac->ams_proj_freq)); 1289 } 1290 } else { 1291 CHKERRQ(PetscViewerASCIIPrintf(viewer," scalar Poisson solver not used (zero-conductivity everywhere) \n")); 1292 } 1293 } 1294 PetscFunctionReturn(0); 1295 } 1296 1297 static PetscErrorCode PCSetFromOptions_HYPRE_ADS(PetscOptionItems *PetscOptionsObject,PC pc) 1298 { 1299 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1300 PetscInt n; 1301 PetscBool flag,flag2,flag3,flag4; 1302 1303 PetscFunctionBegin; 1304 CHKERRQ(PetscOptionsHead(PetscOptionsObject,"HYPRE ADS Options")); 1305 CHKERRQ(PetscOptionsInt("-pc_hypre_ads_print_level","Debugging output level for ADS","None",jac->as_print,&jac->as_print,&flag)); 1306 if (flag) PetscStackCallStandard(HYPRE_ADSSetPrintLevel,jac->hsolver,jac->as_print); 1307 CHKERRQ(PetscOptionsInt("-pc_hypre_ads_max_iter","Maximum number of ADS multigrid iterations within PCApply","None",jac->as_max_iter,&jac->as_max_iter,&flag)); 1308 if (flag) PetscStackCallStandard(HYPRE_ADSSetMaxIter,jac->hsolver,jac->as_max_iter); 1309 CHKERRQ(PetscOptionsInt("-pc_hypre_ads_cycle_type","Cycle type for ADS multigrid","None",jac->ads_cycle_type,&jac->ads_cycle_type,&flag)); 1310 if (flag) PetscStackCallStandard(HYPRE_ADSSetCycleType,jac->hsolver,jac->ads_cycle_type); 1311 CHKERRQ(PetscOptionsReal("-pc_hypre_ads_tol","Error tolerance for ADS multigrid","None",jac->as_tol,&jac->as_tol,&flag)); 1312 if (flag) PetscStackCallStandard(HYPRE_ADSSetTol,jac->hsolver,jac->as_tol); 1313 CHKERRQ(PetscOptionsInt("-pc_hypre_ads_relax_type","Relaxation type for ADS smoother","None",jac->as_relax_type,&jac->as_relax_type,&flag)); 1314 CHKERRQ(PetscOptionsInt("-pc_hypre_ads_relax_times","Number of relaxation steps for ADS smoother","None",jac->as_relax_times,&jac->as_relax_times,&flag2)); 1315 CHKERRQ(PetscOptionsReal("-pc_hypre_ads_relax_weight","Relaxation weight for ADS smoother","None",jac->as_relax_weight,&jac->as_relax_weight,&flag3)); 1316 CHKERRQ(PetscOptionsReal("-pc_hypre_ads_omega","SSOR coefficient for ADS smoother","None",jac->as_omega,&jac->as_omega,&flag4)); 1317 if (flag || flag2 || flag3 || flag4) { 1318 PetscStackCallStandard(HYPRE_ADSSetSmoothingOptions,jac->hsolver,jac->as_relax_type, 1319 jac->as_relax_times, 1320 jac->as_relax_weight, 1321 jac->as_omega); 1322 } 1323 CHKERRQ(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)); 1324 n = 5; 1325 CHKERRQ(PetscOptionsIntArray("-pc_hypre_ads_ams_options","AMG options for AMS solver inside ADS","None",jac->as_amg_alpha_opts,&n,&flag2)); 1326 CHKERRQ(PetscOptionsInt("-pc_hypre_ads_ams_cycle_type","Cycle type for AMS solver inside ADS","None",jac->ams_cycle_type,&jac->ams_cycle_type,&flag3)); 1327 if (flag || flag2 || flag3) { 1328 PetscStackCallStandard(HYPRE_ADSSetAMSOptions,jac->hsolver,jac->ams_cycle_type, /* AMS cycle type */ 1329 jac->as_amg_alpha_opts[0], /* AMG coarsen type */ 1330 jac->as_amg_alpha_opts[1], /* AMG agg_levels */ 1331 jac->as_amg_alpha_opts[2], /* AMG relax_type */ 1332 jac->as_amg_alpha_theta, 1333 jac->as_amg_alpha_opts[3], /* AMG interp_type */ 1334 jac->as_amg_alpha_opts[4]); /* AMG Pmax */ 1335 } 1336 CHKERRQ(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)); 1337 n = 5; 1338 CHKERRQ(PetscOptionsIntArray("-pc_hypre_ads_amg_options","AMG options for vector AMG solver inside ADS","None",jac->as_amg_beta_opts,&n,&flag2)); 1339 if (flag || flag2) { 1340 PetscStackCallStandard(HYPRE_ADSSetAMGOptions,jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */ 1341 jac->as_amg_beta_opts[1], /* AMG agg_levels */ 1342 jac->as_amg_beta_opts[2], /* AMG relax_type */ 1343 jac->as_amg_beta_theta, 1344 jac->as_amg_beta_opts[3], /* AMG interp_type */ 1345 jac->as_amg_beta_opts[4]); /* AMG Pmax */ 1346 } 1347 CHKERRQ(PetscOptionsTail()); 1348 PetscFunctionReturn(0); 1349 } 1350 1351 static PetscErrorCode PCView_HYPRE_ADS(PC pc,PetscViewer viewer) 1352 { 1353 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1354 PetscBool iascii; 1355 1356 PetscFunctionBegin; 1357 CHKERRQ(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii)); 1358 if (iascii) { 1359 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE ADS preconditioning\n")); 1360 CHKERRQ(PetscViewerASCIIPrintf(viewer," subspace iterations per application %d\n",jac->as_max_iter)); 1361 CHKERRQ(PetscViewerASCIIPrintf(viewer," subspace cycle type %d\n",jac->ads_cycle_type)); 1362 CHKERRQ(PetscViewerASCIIPrintf(viewer," subspace iteration tolerance %g\n",jac->as_tol)); 1363 CHKERRQ(PetscViewerASCIIPrintf(viewer," smoother type %d\n",jac->as_relax_type)); 1364 CHKERRQ(PetscViewerASCIIPrintf(viewer," number of smoothing steps %d\n",jac->as_relax_times)); 1365 CHKERRQ(PetscViewerASCIIPrintf(viewer," smoother weight %g\n",jac->as_relax_weight)); 1366 CHKERRQ(PetscViewerASCIIPrintf(viewer," smoother omega %g\n",jac->as_omega)); 1367 CHKERRQ(PetscViewerASCIIPrintf(viewer," AMS solver using boomerAMG\n")); 1368 CHKERRQ(PetscViewerASCIIPrintf(viewer," subspace cycle type %d\n",jac->ams_cycle_type)); 1369 CHKERRQ(PetscViewerASCIIPrintf(viewer," coarsening type %d\n",jac->as_amg_alpha_opts[0])); 1370 CHKERRQ(PetscViewerASCIIPrintf(viewer," levels of aggressive coarsening %d\n",jac->as_amg_alpha_opts[1])); 1371 CHKERRQ(PetscViewerASCIIPrintf(viewer," relaxation type %d\n",jac->as_amg_alpha_opts[2])); 1372 CHKERRQ(PetscViewerASCIIPrintf(viewer," interpolation type %d\n",jac->as_amg_alpha_opts[3])); 1373 CHKERRQ(PetscViewerASCIIPrintf(viewer," max nonzero elements in interpolation rows %d\n",jac->as_amg_alpha_opts[4])); 1374 CHKERRQ(PetscViewerASCIIPrintf(viewer," strength threshold %g\n",jac->as_amg_alpha_theta)); 1375 CHKERRQ(PetscViewerASCIIPrintf(viewer," vector Poisson solver using boomerAMG\n")); 1376 CHKERRQ(PetscViewerASCIIPrintf(viewer," coarsening type %d\n",jac->as_amg_beta_opts[0])); 1377 CHKERRQ(PetscViewerASCIIPrintf(viewer," levels of aggressive coarsening %d\n",jac->as_amg_beta_opts[1])); 1378 CHKERRQ(PetscViewerASCIIPrintf(viewer," relaxation type %d\n",jac->as_amg_beta_opts[2])); 1379 CHKERRQ(PetscViewerASCIIPrintf(viewer," interpolation type %d\n",jac->as_amg_beta_opts[3])); 1380 CHKERRQ(PetscViewerASCIIPrintf(viewer," max nonzero elements in interpolation rows %d\n",jac->as_amg_beta_opts[4])); 1381 CHKERRQ(PetscViewerASCIIPrintf(viewer," strength threshold %g\n",jac->as_amg_beta_theta)); 1382 } 1383 PetscFunctionReturn(0); 1384 } 1385 1386 static PetscErrorCode PCHYPRESetDiscreteGradient_HYPRE(PC pc, Mat G) 1387 { 1388 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1389 PetscBool ishypre; 1390 1391 PetscFunctionBegin; 1392 CHKERRQ(PetscObjectTypeCompare((PetscObject)G,MATHYPRE,&ishypre)); 1393 if (ishypre) { 1394 CHKERRQ(PetscObjectReference((PetscObject)G)); 1395 CHKERRQ(MatDestroy(&jac->G)); 1396 jac->G = G; 1397 } else { 1398 CHKERRQ(MatDestroy(&jac->G)); 1399 CHKERRQ(MatConvert(G,MATHYPRE,MAT_INITIAL_MATRIX,&jac->G)); 1400 } 1401 PetscFunctionReturn(0); 1402 } 1403 1404 /*@ 1405 PCHYPRESetDiscreteGradient - Set discrete gradient matrix 1406 1407 Collective on PC 1408 1409 Input Parameters: 1410 + pc - the preconditioning context 1411 - G - the discrete gradient 1412 1413 Level: intermediate 1414 1415 Notes: 1416 G should have as many rows as the number of edges and as many columns as the number of vertices in the mesh 1417 1418 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 1419 1420 .seealso: PCHYPRESetDiscreteCurl() 1421 @*/ 1422 PetscErrorCode PCHYPRESetDiscreteGradient(PC pc, Mat G) 1423 { 1424 PetscFunctionBegin; 1425 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1426 PetscValidHeaderSpecific(G,MAT_CLASSID,2); 1427 PetscCheckSameComm(pc,1,G,2); 1428 CHKERRQ(PetscTryMethod(pc,"PCHYPRESetDiscreteGradient_C",(PC,Mat),(pc,G))); 1429 PetscFunctionReturn(0); 1430 } 1431 1432 static PetscErrorCode PCHYPRESetDiscreteCurl_HYPRE(PC pc, Mat C) 1433 { 1434 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1435 PetscBool ishypre; 1436 1437 PetscFunctionBegin; 1438 CHKERRQ(PetscObjectTypeCompare((PetscObject)C,MATHYPRE,&ishypre)); 1439 if (ishypre) { 1440 CHKERRQ(PetscObjectReference((PetscObject)C)); 1441 CHKERRQ(MatDestroy(&jac->C)); 1442 jac->C = C; 1443 } else { 1444 CHKERRQ(MatDestroy(&jac->C)); 1445 CHKERRQ(MatConvert(C,MATHYPRE,MAT_INITIAL_MATRIX,&jac->C)); 1446 } 1447 PetscFunctionReturn(0); 1448 } 1449 1450 /*@ 1451 PCHYPRESetDiscreteCurl - Set discrete curl matrix 1452 1453 Collective on PC 1454 1455 Input Parameters: 1456 + pc - the preconditioning context 1457 - C - the discrete curl 1458 1459 Level: intermediate 1460 1461 Notes: 1462 C should have as many rows as the number of faces and as many columns as the number of edges in the mesh 1463 1464 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 1465 1466 .seealso: PCHYPRESetDiscreteGradient() 1467 @*/ 1468 PetscErrorCode PCHYPRESetDiscreteCurl(PC pc, Mat C) 1469 { 1470 PetscFunctionBegin; 1471 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1472 PetscValidHeaderSpecific(C,MAT_CLASSID,2); 1473 PetscCheckSameComm(pc,1,C,2); 1474 CHKERRQ(PetscTryMethod(pc,"PCHYPRESetDiscreteCurl_C",(PC,Mat),(pc,C))); 1475 PetscFunctionReturn(0); 1476 } 1477 1478 static PetscErrorCode PCHYPRESetInterpolations_HYPRE(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[]) 1479 { 1480 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1481 PetscBool ishypre; 1482 PetscInt i; 1483 PetscFunctionBegin; 1484 1485 CHKERRQ(MatDestroy(&jac->RT_PiFull)); 1486 CHKERRQ(MatDestroy(&jac->ND_PiFull)); 1487 for (i=0;i<3;++i) { 1488 CHKERRQ(MatDestroy(&jac->RT_Pi[i])); 1489 CHKERRQ(MatDestroy(&jac->ND_Pi[i])); 1490 } 1491 1492 jac->dim = dim; 1493 if (RT_PiFull) { 1494 CHKERRQ(PetscObjectTypeCompare((PetscObject)RT_PiFull,MATHYPRE,&ishypre)); 1495 if (ishypre) { 1496 CHKERRQ(PetscObjectReference((PetscObject)RT_PiFull)); 1497 jac->RT_PiFull = RT_PiFull; 1498 } else { 1499 CHKERRQ(MatConvert(RT_PiFull,MATHYPRE,MAT_INITIAL_MATRIX,&jac->RT_PiFull)); 1500 } 1501 } 1502 if (RT_Pi) { 1503 for (i=0;i<dim;++i) { 1504 if (RT_Pi[i]) { 1505 CHKERRQ(PetscObjectTypeCompare((PetscObject)RT_Pi[i],MATHYPRE,&ishypre)); 1506 if (ishypre) { 1507 CHKERRQ(PetscObjectReference((PetscObject)RT_Pi[i])); 1508 jac->RT_Pi[i] = RT_Pi[i]; 1509 } else { 1510 CHKERRQ(MatConvert(RT_Pi[i],MATHYPRE,MAT_INITIAL_MATRIX,&jac->RT_Pi[i])); 1511 } 1512 } 1513 } 1514 } 1515 if (ND_PiFull) { 1516 CHKERRQ(PetscObjectTypeCompare((PetscObject)ND_PiFull,MATHYPRE,&ishypre)); 1517 if (ishypre) { 1518 CHKERRQ(PetscObjectReference((PetscObject)ND_PiFull)); 1519 jac->ND_PiFull = ND_PiFull; 1520 } else { 1521 CHKERRQ(MatConvert(ND_PiFull,MATHYPRE,MAT_INITIAL_MATRIX,&jac->ND_PiFull)); 1522 } 1523 } 1524 if (ND_Pi) { 1525 for (i=0;i<dim;++i) { 1526 if (ND_Pi[i]) { 1527 CHKERRQ(PetscObjectTypeCompare((PetscObject)ND_Pi[i],MATHYPRE,&ishypre)); 1528 if (ishypre) { 1529 CHKERRQ(PetscObjectReference((PetscObject)ND_Pi[i])); 1530 jac->ND_Pi[i] = ND_Pi[i]; 1531 } else { 1532 CHKERRQ(MatConvert(ND_Pi[i],MATHYPRE,MAT_INITIAL_MATRIX,&jac->ND_Pi[i])); 1533 } 1534 } 1535 } 1536 } 1537 1538 PetscFunctionReturn(0); 1539 } 1540 1541 /*@ 1542 PCHYPRESetInterpolations - Set interpolation matrices for AMS/ADS preconditioner 1543 1544 Collective on PC 1545 1546 Input Parameters: 1547 + pc - the preconditioning context 1548 - dim - the dimension of the problem, only used in AMS 1549 - RT_PiFull - Raviart-Thomas interpolation matrix 1550 - RT_Pi - x/y/z component of Raviart-Thomas interpolation matrix 1551 - ND_PiFull - Nedelec interpolation matrix 1552 - ND_Pi - x/y/z component of Nedelec interpolation matrix 1553 1554 Notes: 1555 For AMS, only Nedelec interpolation matrices are needed, the Raviart-Thomas interpolation matrices can be set to NULL. 1556 1557 For ADS, both type of interpolation matrices are needed. 1558 1559 Level: intermediate 1560 1561 @*/ 1562 PetscErrorCode PCHYPRESetInterpolations(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[]) 1563 { 1564 PetscInt i; 1565 1566 PetscFunctionBegin; 1567 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1568 if (RT_PiFull) { 1569 PetscValidHeaderSpecific(RT_PiFull,MAT_CLASSID,3); 1570 PetscCheckSameComm(pc,1,RT_PiFull,3); 1571 } 1572 if (RT_Pi) { 1573 PetscValidPointer(RT_Pi,4); 1574 for (i=0;i<dim;++i) { 1575 if (RT_Pi[i]) { 1576 PetscValidHeaderSpecific(RT_Pi[i],MAT_CLASSID,4); 1577 PetscCheckSameComm(pc,1,RT_Pi[i],4); 1578 } 1579 } 1580 } 1581 if (ND_PiFull) { 1582 PetscValidHeaderSpecific(ND_PiFull,MAT_CLASSID,5); 1583 PetscCheckSameComm(pc,1,ND_PiFull,5); 1584 } 1585 if (ND_Pi) { 1586 PetscValidPointer(ND_Pi,6); 1587 for (i=0;i<dim;++i) { 1588 if (ND_Pi[i]) { 1589 PetscValidHeaderSpecific(ND_Pi[i],MAT_CLASSID,6); 1590 PetscCheckSameComm(pc,1,ND_Pi[i],6); 1591 } 1592 } 1593 } 1594 CHKERRQ(PetscTryMethod(pc,"PCHYPRESetInterpolations_C",(PC,PetscInt,Mat,Mat[],Mat,Mat[]),(pc,dim,RT_PiFull,RT_Pi,ND_PiFull,ND_Pi))); 1595 PetscFunctionReturn(0); 1596 } 1597 1598 static PetscErrorCode PCHYPRESetPoissonMatrix_HYPRE(PC pc, Mat A, PetscBool isalpha) 1599 { 1600 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1601 PetscBool ishypre; 1602 1603 PetscFunctionBegin; 1604 CHKERRQ(PetscObjectTypeCompare((PetscObject)A,MATHYPRE,&ishypre)); 1605 if (ishypre) { 1606 if (isalpha) { 1607 CHKERRQ(PetscObjectReference((PetscObject)A)); 1608 CHKERRQ(MatDestroy(&jac->alpha_Poisson)); 1609 jac->alpha_Poisson = A; 1610 } else { 1611 if (A) { 1612 CHKERRQ(PetscObjectReference((PetscObject)A)); 1613 } else { 1614 jac->ams_beta_is_zero = PETSC_TRUE; 1615 } 1616 CHKERRQ(MatDestroy(&jac->beta_Poisson)); 1617 jac->beta_Poisson = A; 1618 } 1619 } else { 1620 if (isalpha) { 1621 CHKERRQ(MatDestroy(&jac->alpha_Poisson)); 1622 CHKERRQ(MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&jac->alpha_Poisson)); 1623 } else { 1624 if (A) { 1625 CHKERRQ(MatDestroy(&jac->beta_Poisson)); 1626 CHKERRQ(MatConvert(A,MATHYPRE,MAT_INITIAL_MATRIX,&jac->beta_Poisson)); 1627 } else { 1628 CHKERRQ(MatDestroy(&jac->beta_Poisson)); 1629 jac->ams_beta_is_zero = PETSC_TRUE; 1630 } 1631 } 1632 } 1633 PetscFunctionReturn(0); 1634 } 1635 1636 /*@ 1637 PCHYPRESetAlphaPoissonMatrix - Set vector Poisson matrix 1638 1639 Collective on PC 1640 1641 Input Parameters: 1642 + pc - the preconditioning context 1643 - A - the matrix 1644 1645 Level: intermediate 1646 1647 Notes: 1648 A should be obtained by discretizing the vector valued Poisson problem with linear finite elements 1649 1650 .seealso: PCHYPRESetDiscreteGradient(), PCHYPRESetDiscreteCurl(), PCHYPRESetBetaPoissonMatrix() 1651 @*/ 1652 PetscErrorCode PCHYPRESetAlphaPoissonMatrix(PC pc, Mat A) 1653 { 1654 PetscFunctionBegin; 1655 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1656 PetscValidHeaderSpecific(A,MAT_CLASSID,2); 1657 PetscCheckSameComm(pc,1,A,2); 1658 CHKERRQ(PetscTryMethod(pc,"PCHYPRESetPoissonMatrix_C",(PC,Mat,PetscBool),(pc,A,PETSC_TRUE))); 1659 PetscFunctionReturn(0); 1660 } 1661 1662 /*@ 1663 PCHYPRESetBetaPoissonMatrix - Set Poisson matrix 1664 1665 Collective on PC 1666 1667 Input Parameters: 1668 + pc - the preconditioning context 1669 - A - the matrix 1670 1671 Level: intermediate 1672 1673 Notes: 1674 A should be obtained by discretizing the Poisson problem with linear finite elements. 1675 Following HYPRE convention, the scalar Poisson solver of AMS can be turned off by passing NULL. 1676 1677 .seealso: PCHYPRESetDiscreteGradient(), PCHYPRESetDiscreteCurl(), PCHYPRESetAlphaPoissonMatrix() 1678 @*/ 1679 PetscErrorCode PCHYPRESetBetaPoissonMatrix(PC pc, Mat A) 1680 { 1681 PetscFunctionBegin; 1682 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1683 if (A) { 1684 PetscValidHeaderSpecific(A,MAT_CLASSID,2); 1685 PetscCheckSameComm(pc,1,A,2); 1686 } 1687 CHKERRQ(PetscTryMethod(pc,"PCHYPRESetPoissonMatrix_C",(PC,Mat,PetscBool),(pc,A,PETSC_FALSE))); 1688 PetscFunctionReturn(0); 1689 } 1690 1691 static PetscErrorCode PCHYPRESetEdgeConstantVectors_HYPRE(PC pc,Vec ozz, Vec zoz, Vec zzo) 1692 { 1693 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1694 1695 PetscFunctionBegin; 1696 /* throw away any vector if already set */ 1697 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->constants[0])); 1698 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->constants[1])); 1699 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->constants[2])); 1700 CHKERRQ(VecHYPRE_IJVectorCreate(ozz->map,&jac->constants[0])); 1701 CHKERRQ(VecHYPRE_IJVectorCopy(ozz,jac->constants[0])); 1702 CHKERRQ(VecHYPRE_IJVectorCreate(zoz->map,&jac->constants[1])); 1703 CHKERRQ(VecHYPRE_IJVectorCopy(zoz,jac->constants[1])); 1704 jac->dim = 2; 1705 if (zzo) { 1706 CHKERRQ(VecHYPRE_IJVectorCreate(zzo->map,&jac->constants[2])); 1707 CHKERRQ(VecHYPRE_IJVectorCopy(zzo,jac->constants[2])); 1708 jac->dim++; 1709 } 1710 PetscFunctionReturn(0); 1711 } 1712 1713 /*@ 1714 PCHYPRESetEdgeConstantVectors - Set the representation of the constant vector fields in the edge element basis 1715 1716 Collective on PC 1717 1718 Input Parameters: 1719 + pc - the preconditioning context 1720 - ozz - vector representing (1,0,0) (or (1,0) in 2D) 1721 - zoz - vector representing (0,1,0) (or (0,1) in 2D) 1722 - zzo - vector representing (0,0,1) (use NULL in 2D) 1723 1724 Level: intermediate 1725 1726 @*/ 1727 PetscErrorCode PCHYPRESetEdgeConstantVectors(PC pc, Vec ozz, Vec zoz, Vec zzo) 1728 { 1729 PetscFunctionBegin; 1730 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1731 PetscValidHeaderSpecific(ozz,VEC_CLASSID,2); 1732 PetscValidHeaderSpecific(zoz,VEC_CLASSID,3); 1733 if (zzo) PetscValidHeaderSpecific(zzo,VEC_CLASSID,4); 1734 PetscCheckSameComm(pc,1,ozz,2); 1735 PetscCheckSameComm(pc,1,zoz,3); 1736 if (zzo) PetscCheckSameComm(pc,1,zzo,4); 1737 CHKERRQ(PetscTryMethod(pc,"PCHYPRESetEdgeConstantVectors_C",(PC,Vec,Vec,Vec),(pc,ozz,zoz,zzo))); 1738 PetscFunctionReturn(0); 1739 } 1740 1741 static PetscErrorCode PCSetCoordinates_HYPRE(PC pc, PetscInt dim, PetscInt nloc, PetscReal *coords) 1742 { 1743 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1744 Vec tv; 1745 PetscInt i; 1746 1747 PetscFunctionBegin; 1748 /* throw away any coordinate vector if already set */ 1749 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->coords[0])); 1750 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->coords[1])); 1751 CHKERRQ(VecHYPRE_IJVectorDestroy(&jac->coords[2])); 1752 jac->dim = dim; 1753 1754 /* compute IJ vector for coordinates */ 1755 CHKERRQ(VecCreate(PetscObjectComm((PetscObject)pc),&tv)); 1756 CHKERRQ(VecSetType(tv,VECSTANDARD)); 1757 CHKERRQ(VecSetSizes(tv,nloc,PETSC_DECIDE)); 1758 for (i=0;i<dim;i++) { 1759 PetscScalar *array; 1760 PetscInt j; 1761 1762 CHKERRQ(VecHYPRE_IJVectorCreate(tv->map,&jac->coords[i])); 1763 CHKERRQ(VecGetArrayWrite(tv,&array)); 1764 for (j=0;j<nloc;j++) array[j] = coords[j*dim+i]; 1765 CHKERRQ(VecRestoreArrayWrite(tv,&array)); 1766 CHKERRQ(VecHYPRE_IJVectorCopy(tv,jac->coords[i])); 1767 } 1768 CHKERRQ(VecDestroy(&tv)); 1769 PetscFunctionReturn(0); 1770 } 1771 1772 /* ---------------------------------------------------------------------------------*/ 1773 1774 static PetscErrorCode PCHYPREGetType_HYPRE(PC pc,const char *name[]) 1775 { 1776 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1777 1778 PetscFunctionBegin; 1779 *name = jac->hypre_type; 1780 PetscFunctionReturn(0); 1781 } 1782 1783 static PetscErrorCode PCHYPRESetType_HYPRE(PC pc,const char name[]) 1784 { 1785 PC_HYPRE *jac = (PC_HYPRE*)pc->data; 1786 PetscBool flag; 1787 1788 PetscFunctionBegin; 1789 if (jac->hypre_type) { 1790 CHKERRQ(PetscStrcmp(jac->hypre_type,name,&flag)); 1791 PetscCheck(flag,PetscObjectComm((PetscObject)pc),PETSC_ERR_ORDER,"Cannot reset the HYPRE preconditioner type once it has been set"); 1792 PetscFunctionReturn(0); 1793 } else { 1794 CHKERRQ(PetscStrallocpy(name, &jac->hypre_type)); 1795 } 1796 1797 jac->maxiter = PETSC_DEFAULT; 1798 jac->tol = PETSC_DEFAULT; 1799 jac->printstatistics = PetscLogPrintInfo; 1800 1801 CHKERRQ(PetscStrcmp("pilut",jac->hypre_type,&flag)); 1802 if (flag) { 1803 CHKERRQ(PetscCommGetComm(PetscObjectComm((PetscObject)pc),&jac->comm_hypre)); 1804 PetscStackCallStandard(HYPRE_ParCSRPilutCreate,jac->comm_hypre,&jac->hsolver); 1805 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Pilut; 1806 pc->ops->view = PCView_HYPRE_Pilut; 1807 jac->destroy = HYPRE_ParCSRPilutDestroy; 1808 jac->setup = HYPRE_ParCSRPilutSetup; 1809 jac->solve = HYPRE_ParCSRPilutSolve; 1810 jac->factorrowsize = PETSC_DEFAULT; 1811 PetscFunctionReturn(0); 1812 } 1813 CHKERRQ(PetscStrcmp("euclid",jac->hypre_type,&flag)); 1814 if (flag) { 1815 #if defined(PETSC_HAVE_64BIT_INDICES) 1816 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Hypre Euclid not support with 64 bit indices"); 1817 #endif 1818 CHKERRQ(PetscCommGetComm(PetscObjectComm((PetscObject)pc),&jac->comm_hypre)); 1819 PetscStackCallStandard(HYPRE_EuclidCreate,jac->comm_hypre,&jac->hsolver); 1820 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Euclid; 1821 pc->ops->view = PCView_HYPRE_Euclid; 1822 jac->destroy = HYPRE_EuclidDestroy; 1823 jac->setup = HYPRE_EuclidSetup; 1824 jac->solve = HYPRE_EuclidSolve; 1825 jac->factorrowsize = PETSC_DEFAULT; 1826 jac->eu_level = PETSC_DEFAULT; /* default */ 1827 PetscFunctionReturn(0); 1828 } 1829 CHKERRQ(PetscStrcmp("parasails",jac->hypre_type,&flag)); 1830 if (flag) { 1831 CHKERRQ(PetscCommGetComm(PetscObjectComm((PetscObject)pc),&jac->comm_hypre)); 1832 PetscStackCallStandard(HYPRE_ParaSailsCreate,jac->comm_hypre,&jac->hsolver); 1833 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ParaSails; 1834 pc->ops->view = PCView_HYPRE_ParaSails; 1835 jac->destroy = HYPRE_ParaSailsDestroy; 1836 jac->setup = HYPRE_ParaSailsSetup; 1837 jac->solve = HYPRE_ParaSailsSolve; 1838 /* initialize */ 1839 jac->nlevels = 1; 1840 jac->threshold = .1; 1841 jac->filter = .1; 1842 jac->loadbal = 0; 1843 if (PetscLogPrintInfo) jac->logging = (int) PETSC_TRUE; 1844 else jac->logging = (int) PETSC_FALSE; 1845 1846 jac->ruse = (int) PETSC_FALSE; 1847 jac->symt = 0; 1848 PetscStackCallStandard(HYPRE_ParaSailsSetParams,jac->hsolver,jac->threshold,jac->nlevels); 1849 PetscStackCallStandard(HYPRE_ParaSailsSetFilter,jac->hsolver,jac->filter); 1850 PetscStackCallStandard(HYPRE_ParaSailsSetLoadbal,jac->hsolver,jac->loadbal); 1851 PetscStackCallStandard(HYPRE_ParaSailsSetLogging,jac->hsolver,jac->logging); 1852 PetscStackCallStandard(HYPRE_ParaSailsSetReuse,jac->hsolver,jac->ruse); 1853 PetscStackCallStandard(HYPRE_ParaSailsSetSym,jac->hsolver,jac->symt); 1854 PetscFunctionReturn(0); 1855 } 1856 CHKERRQ(PetscStrcmp("boomeramg",jac->hypre_type,&flag)); 1857 if (flag) { 1858 PetscStackCallStandard(HYPRE_BoomerAMGCreate,&jac->hsolver); 1859 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_BoomerAMG; 1860 pc->ops->view = PCView_HYPRE_BoomerAMG; 1861 pc->ops->applytranspose = PCApplyTranspose_HYPRE_BoomerAMG; 1862 pc->ops->applyrichardson = PCApplyRichardson_HYPRE_BoomerAMG; 1863 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCGetInterpolations_C",PCGetInterpolations_BoomerAMG)); 1864 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCGetCoarseOperators_C",PCGetCoarseOperators_BoomerAMG)); 1865 jac->destroy = HYPRE_BoomerAMGDestroy; 1866 jac->setup = HYPRE_BoomerAMGSetup; 1867 jac->solve = HYPRE_BoomerAMGSolve; 1868 jac->applyrichardson = PETSC_FALSE; 1869 /* these defaults match the hypre defaults */ 1870 jac->cycletype = 1; 1871 jac->maxlevels = 25; 1872 jac->maxiter = 1; 1873 jac->tol = 0.0; /* tolerance of zero indicates use as preconditioner (suppresses convergence errors) */ 1874 jac->truncfactor = 0.0; 1875 jac->strongthreshold = .25; 1876 jac->maxrowsum = .9; 1877 jac->coarsentype = 6; 1878 jac->measuretype = 0; 1879 jac->gridsweeps[0] = jac->gridsweeps[1] = jac->gridsweeps[2] = 1; 1880 jac->smoothtype = -1; /* Not set by default */ 1881 jac->smoothnumlevels = 25; 1882 jac->eu_level = 0; 1883 jac->eu_droptolerance = 0; 1884 jac->eu_bj = 0; 1885 jac->relaxtype[0] = jac->relaxtype[1] = 6; /* Defaults to SYMMETRIC since in PETSc we are using a PC - most likely with CG */ 1886 jac->relaxtype[2] = 9; /*G.E. */ 1887 jac->relaxweight = 1.0; 1888 jac->outerrelaxweight = 1.0; 1889 jac->relaxorder = 1; 1890 jac->interptype = 0; 1891 jac->Rtype = 0; 1892 jac->Rstrongthreshold = 0.25; 1893 jac->Rfilterthreshold = 0.0; 1894 jac->Adroptype = -1; 1895 jac->Adroptol = 0.0; 1896 jac->agg_nl = 0; 1897 jac->agg_interptype = 4; 1898 jac->pmax = 0; 1899 jac->truncfactor = 0.0; 1900 jac->agg_num_paths = 1; 1901 jac->maxc = 9; 1902 jac->minc = 1; 1903 jac->nodal_coarsening = 0; 1904 jac->nodal_coarsening_diag = 0; 1905 jac->vec_interp_variant = 0; 1906 jac->vec_interp_qmax = 0; 1907 jac->vec_interp_smooth = PETSC_FALSE; 1908 jac->interp_refine = 0; 1909 jac->nodal_relax = PETSC_FALSE; 1910 jac->nodal_relax_levels = 1; 1911 jac->rap2 = 0; 1912 1913 /* GPU defaults 1914 from https://hypre.readthedocs.io/en/latest/solvers-boomeramg.html#gpu-supported-options 1915 and /src/parcsr_ls/par_amg.c */ 1916 #if defined(PETSC_HAVE_HYPRE_DEVICE) 1917 jac->keeptranspose = PETSC_TRUE; 1918 jac->mod_rap2 = 1; 1919 jac->coarsentype = 8; 1920 jac->relaxorder = 0; 1921 jac->interptype = 6; 1922 jac->relaxtype[0] = 18; 1923 jac->relaxtype[1] = 18; 1924 jac->agg_interptype = 7; 1925 #else 1926 jac->keeptranspose = PETSC_FALSE; 1927 jac->mod_rap2 = 0; 1928 #endif 1929 PetscStackCallStandard(HYPRE_BoomerAMGSetCycleType,jac->hsolver,jac->cycletype); 1930 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxLevels,jac->hsolver,jac->maxlevels); 1931 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxIter,jac->hsolver,jac->maxiter); 1932 PetscStackCallStandard(HYPRE_BoomerAMGSetTol,jac->hsolver,jac->tol); 1933 PetscStackCallStandard(HYPRE_BoomerAMGSetTruncFactor,jac->hsolver,jac->truncfactor); 1934 PetscStackCallStandard(HYPRE_BoomerAMGSetStrongThreshold,jac->hsolver,jac->strongthreshold); 1935 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxRowSum,jac->hsolver,jac->maxrowsum); 1936 PetscStackCallStandard(HYPRE_BoomerAMGSetCoarsenType,jac->hsolver,jac->coarsentype); 1937 PetscStackCallStandard(HYPRE_BoomerAMGSetMeasureType,jac->hsolver,jac->measuretype); 1938 PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxOrder,jac->hsolver, jac->relaxorder); 1939 PetscStackCallStandard(HYPRE_BoomerAMGSetInterpType,jac->hsolver,jac->interptype); 1940 PetscStackCallStandard(HYPRE_BoomerAMGSetAggNumLevels,jac->hsolver,jac->agg_nl); 1941 PetscStackCallStandard(HYPRE_BoomerAMGSetAggInterpType,jac->hsolver,jac->agg_interptype); 1942 PetscStackCallStandard(HYPRE_BoomerAMGSetPMaxElmts,jac->hsolver,jac->pmax); 1943 PetscStackCallStandard(HYPRE_BoomerAMGSetNumPaths,jac->hsolver,jac->agg_num_paths); 1944 PetscStackCallStandard(HYPRE_BoomerAMGSetRelaxType,jac->hsolver, jac->relaxtype[0]); /* defaults coarse to 9 */ 1945 PetscStackCallStandard(HYPRE_BoomerAMGSetNumSweeps,jac->hsolver, jac->gridsweeps[0]); /* defaults coarse to 1 */ 1946 PetscStackCallStandard(HYPRE_BoomerAMGSetMaxCoarseSize,jac->hsolver, jac->maxc); 1947 PetscStackCallStandard(HYPRE_BoomerAMGSetMinCoarseSize,jac->hsolver, jac->minc); 1948 /* GPU */ 1949 #if PETSC_PKG_HYPRE_VERSION_GE(2,18,0) 1950 PetscStackCallStandard(HYPRE_BoomerAMGSetKeepTranspose,jac->hsolver,jac->keeptranspose ? 1 : 0); 1951 PetscStackCallStandard(HYPRE_BoomerAMGSetRAP2,jac->hsolver, jac->rap2); 1952 PetscStackCallStandard(HYPRE_BoomerAMGSetModuleRAP2,jac->hsolver, jac->mod_rap2); 1953 #endif 1954 1955 /* AIR */ 1956 #if PETSC_PKG_HYPRE_VERSION_GE(2,18,0) 1957 PetscStackCallStandard(HYPRE_BoomerAMGSetRestriction,jac->hsolver,jac->Rtype); 1958 PetscStackCallStandard(HYPRE_BoomerAMGSetStrongThresholdR,jac->hsolver,jac->Rstrongthreshold); 1959 PetscStackCallStandard(HYPRE_BoomerAMGSetFilterThresholdR,jac->hsolver,jac->Rfilterthreshold); 1960 PetscStackCallStandard(HYPRE_BoomerAMGSetADropTol,jac->hsolver,jac->Adroptol); 1961 PetscStackCallStandard(HYPRE_BoomerAMGSetADropType,jac->hsolver,jac->Adroptype); 1962 #endif 1963 PetscFunctionReturn(0); 1964 } 1965 CHKERRQ(PetscStrcmp("ams",jac->hypre_type,&flag)); 1966 if (flag) { 1967 CHKERRQ(HYPRE_AMSCreate(&jac->hsolver)); 1968 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_AMS; 1969 pc->ops->view = PCView_HYPRE_AMS; 1970 jac->destroy = HYPRE_AMSDestroy; 1971 jac->setup = HYPRE_AMSSetup; 1972 jac->solve = HYPRE_AMSSolve; 1973 jac->coords[0] = NULL; 1974 jac->coords[1] = NULL; 1975 jac->coords[2] = NULL; 1976 /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE AMS */ 1977 jac->as_print = 0; 1978 jac->as_max_iter = 1; /* used as a preconditioner */ 1979 jac->as_tol = 0.; /* used as a preconditioner */ 1980 jac->ams_cycle_type = 13; 1981 /* Smoothing options */ 1982 jac->as_relax_type = 2; 1983 jac->as_relax_times = 1; 1984 jac->as_relax_weight = 1.0; 1985 jac->as_omega = 1.0; 1986 /* Vector valued Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */ 1987 jac->as_amg_alpha_opts[0] = 10; 1988 jac->as_amg_alpha_opts[1] = 1; 1989 jac->as_amg_alpha_opts[2] = 6; 1990 jac->as_amg_alpha_opts[3] = 6; 1991 jac->as_amg_alpha_opts[4] = 4; 1992 jac->as_amg_alpha_theta = 0.25; 1993 /* Scalar Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */ 1994 jac->as_amg_beta_opts[0] = 10; 1995 jac->as_amg_beta_opts[1] = 1; 1996 jac->as_amg_beta_opts[2] = 6; 1997 jac->as_amg_beta_opts[3] = 6; 1998 jac->as_amg_beta_opts[4] = 4; 1999 jac->as_amg_beta_theta = 0.25; 2000 PetscStackCallStandard(HYPRE_AMSSetPrintLevel,jac->hsolver,jac->as_print); 2001 PetscStackCallStandard(HYPRE_AMSSetMaxIter,jac->hsolver,jac->as_max_iter); 2002 PetscStackCallStandard(HYPRE_AMSSetCycleType,jac->hsolver,jac->ams_cycle_type); 2003 PetscStackCallStandard(HYPRE_AMSSetTol,jac->hsolver,jac->as_tol); 2004 PetscStackCallStandard(HYPRE_AMSSetSmoothingOptions,jac->hsolver,jac->as_relax_type, 2005 jac->as_relax_times, 2006 jac->as_relax_weight, 2007 jac->as_omega); 2008 PetscStackCallStandard(HYPRE_AMSSetAlphaAMGOptions,jac->hsolver,jac->as_amg_alpha_opts[0], /* AMG coarsen type */ 2009 jac->as_amg_alpha_opts[1], /* AMG agg_levels */ 2010 jac->as_amg_alpha_opts[2], /* AMG relax_type */ 2011 jac->as_amg_alpha_theta, 2012 jac->as_amg_alpha_opts[3], /* AMG interp_type */ 2013 jac->as_amg_alpha_opts[4]); /* AMG Pmax */ 2014 PetscStackCallStandard(HYPRE_AMSSetBetaAMGOptions,jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */ 2015 jac->as_amg_beta_opts[1], /* AMG agg_levels */ 2016 jac->as_amg_beta_opts[2], /* AMG relax_type */ 2017 jac->as_amg_beta_theta, 2018 jac->as_amg_beta_opts[3], /* AMG interp_type */ 2019 jac->as_amg_beta_opts[4]); /* AMG Pmax */ 2020 /* Zero conductivity */ 2021 jac->ams_beta_is_zero = PETSC_FALSE; 2022 jac->ams_beta_is_zero_part = PETSC_FALSE; 2023 PetscFunctionReturn(0); 2024 } 2025 CHKERRQ(PetscStrcmp("ads",jac->hypre_type,&flag)); 2026 if (flag) { 2027 CHKERRQ(HYPRE_ADSCreate(&jac->hsolver)); 2028 pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ADS; 2029 pc->ops->view = PCView_HYPRE_ADS; 2030 jac->destroy = HYPRE_ADSDestroy; 2031 jac->setup = HYPRE_ADSSetup; 2032 jac->solve = HYPRE_ADSSolve; 2033 jac->coords[0] = NULL; 2034 jac->coords[1] = NULL; 2035 jac->coords[2] = NULL; 2036 /* solver parameters: these are borrowed from mfem package, and they are not the default values from HYPRE ADS */ 2037 jac->as_print = 0; 2038 jac->as_max_iter = 1; /* used as a preconditioner */ 2039 jac->as_tol = 0.; /* used as a preconditioner */ 2040 jac->ads_cycle_type = 13; 2041 /* Smoothing options */ 2042 jac->as_relax_type = 2; 2043 jac->as_relax_times = 1; 2044 jac->as_relax_weight = 1.0; 2045 jac->as_omega = 1.0; 2046 /* AMS solver parameters: cycle_type, coarsen type, agg_levels, relax_type, interp_type, Pmax */ 2047 jac->ams_cycle_type = 14; 2048 jac->as_amg_alpha_opts[0] = 10; 2049 jac->as_amg_alpha_opts[1] = 1; 2050 jac->as_amg_alpha_opts[2] = 6; 2051 jac->as_amg_alpha_opts[3] = 6; 2052 jac->as_amg_alpha_opts[4] = 4; 2053 jac->as_amg_alpha_theta = 0.25; 2054 /* Vector Poisson AMG solver parameters: coarsen type, agg_levels, relax_type, interp_type, Pmax */ 2055 jac->as_amg_beta_opts[0] = 10; 2056 jac->as_amg_beta_opts[1] = 1; 2057 jac->as_amg_beta_opts[2] = 6; 2058 jac->as_amg_beta_opts[3] = 6; 2059 jac->as_amg_beta_opts[4] = 4; 2060 jac->as_amg_beta_theta = 0.25; 2061 PetscStackCallStandard(HYPRE_ADSSetPrintLevel,jac->hsolver,jac->as_print); 2062 PetscStackCallStandard(HYPRE_ADSSetMaxIter,jac->hsolver,jac->as_max_iter); 2063 PetscStackCallStandard(HYPRE_ADSSetCycleType,jac->hsolver,jac->ams_cycle_type); 2064 PetscStackCallStandard(HYPRE_ADSSetTol,jac->hsolver,jac->as_tol); 2065 PetscStackCallStandard(HYPRE_ADSSetSmoothingOptions,jac->hsolver,jac->as_relax_type, 2066 jac->as_relax_times, 2067 jac->as_relax_weight, 2068 jac->as_omega); 2069 PetscStackCallStandard(HYPRE_ADSSetAMSOptions,jac->hsolver,jac->ams_cycle_type, /* AMG coarsen type */ 2070 jac->as_amg_alpha_opts[0], /* AMG coarsen type */ 2071 jac->as_amg_alpha_opts[1], /* AMG agg_levels */ 2072 jac->as_amg_alpha_opts[2], /* AMG relax_type */ 2073 jac->as_amg_alpha_theta, 2074 jac->as_amg_alpha_opts[3], /* AMG interp_type */ 2075 jac->as_amg_alpha_opts[4]); /* AMG Pmax */ 2076 PetscStackCallStandard(HYPRE_ADSSetAMGOptions,jac->hsolver,jac->as_amg_beta_opts[0], /* AMG coarsen type */ 2077 jac->as_amg_beta_opts[1], /* AMG agg_levels */ 2078 jac->as_amg_beta_opts[2], /* AMG relax_type */ 2079 jac->as_amg_beta_theta, 2080 jac->as_amg_beta_opts[3], /* AMG interp_type */ 2081 jac->as_amg_beta_opts[4]); /* AMG Pmax */ 2082 PetscFunctionReturn(0); 2083 } 2084 CHKERRQ(PetscFree(jac->hypre_type)); 2085 2086 jac->hypre_type = NULL; 2087 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown HYPRE preconditioner %s; Choices are euclid, pilut, parasails, boomeramg, ams",name); 2088 } 2089 2090 /* 2091 It only gets here if the HYPRE type has not been set before the call to 2092 ...SetFromOptions() which actually is most of the time 2093 */ 2094 PetscErrorCode PCSetFromOptions_HYPRE(PetscOptionItems *PetscOptionsObject,PC pc) 2095 { 2096 PetscInt indx; 2097 const char *type[] = {"euclid","pilut","parasails","boomeramg","ams","ads"}; 2098 PetscBool flg; 2099 2100 PetscFunctionBegin; 2101 CHKERRQ(PetscOptionsHead(PetscOptionsObject,"HYPRE preconditioner options")); 2102 CHKERRQ(PetscOptionsEList("-pc_hypre_type","HYPRE preconditioner type","PCHYPRESetType",type,ALEN(type),"boomeramg",&indx,&flg)); 2103 if (flg) { 2104 CHKERRQ(PCHYPRESetType_HYPRE(pc,type[indx])); 2105 } else { 2106 CHKERRQ(PCHYPRESetType_HYPRE(pc,"boomeramg")); 2107 } 2108 if (pc->ops->setfromoptions) { 2109 CHKERRQ(pc->ops->setfromoptions(PetscOptionsObject,pc)); 2110 } 2111 CHKERRQ(PetscOptionsTail()); 2112 PetscFunctionReturn(0); 2113 } 2114 2115 /*@C 2116 PCHYPRESetType - Sets which hypre preconditioner you wish to use 2117 2118 Input Parameters: 2119 + pc - the preconditioner context 2120 - name - either euclid, pilut, parasails, boomeramg, ams, ads 2121 2122 Options Database Keys: 2123 -pc_hypre_type - One of euclid, pilut, parasails, boomeramg, ams, ads 2124 2125 Level: intermediate 2126 2127 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 2128 PCHYPRE 2129 2130 @*/ 2131 PetscErrorCode PCHYPRESetType(PC pc,const char name[]) 2132 { 2133 PetscFunctionBegin; 2134 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 2135 PetscValidCharPointer(name,2); 2136 CHKERRQ(PetscTryMethod(pc,"PCHYPRESetType_C",(PC,const char[]),(pc,name))); 2137 PetscFunctionReturn(0); 2138 } 2139 2140 /*@C 2141 PCHYPREGetType - Gets which hypre preconditioner you are using 2142 2143 Input Parameter: 2144 . pc - the preconditioner context 2145 2146 Output Parameter: 2147 . name - either euclid, pilut, parasails, boomeramg, ams, ads 2148 2149 Level: intermediate 2150 2151 .seealso: PCCreate(), PCHYPRESetType(), PCType (for list of available types), PC, 2152 PCHYPRE 2153 2154 @*/ 2155 PetscErrorCode PCHYPREGetType(PC pc,const char *name[]) 2156 { 2157 PetscFunctionBegin; 2158 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 2159 PetscValidPointer(name,2); 2160 CHKERRQ(PetscTryMethod(pc,"PCHYPREGetType_C",(PC,const char*[]),(pc,name))); 2161 PetscFunctionReturn(0); 2162 } 2163 2164 /*@C 2165 PCMGGalerkinSetMatProductAlgorithm - Set type of SpGEMM for hypre to use 2166 2167 Logically Collective on PC 2168 2169 Input Parameters: 2170 + pc - the hypre context 2171 - type - one of 'cusparse', 'hypre' 2172 2173 Options Database Key: 2174 . -pc_mg_galerkin_mat_product_algorithm <cusparse,hypre> - Type of SpGEMM to use in hypre 2175 2176 Level: intermediate 2177 2178 .seealso: PCMGGalerkinGetMatProductAlgorithm() 2179 2180 @*/ 2181 PetscErrorCode PCMGGalerkinSetMatProductAlgorithm(PC pc,const char name[]) 2182 { 2183 PetscFunctionBegin; 2184 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 2185 CHKERRQ(PetscTryMethod(pc,"PCMGGalerkinSetMatProductAlgorithm_C",(PC,const char[]),(pc,name))); 2186 PetscFunctionReturn(0); 2187 } 2188 2189 /*@C 2190 PCMGGalerkinGetMatProductAlgorithm - Get type of SpGEMM for hypre 2191 2192 Not Collective 2193 2194 Input Parameter: 2195 . pc - the multigrid context 2196 2197 Output Parameter: 2198 . name - one of 'cusparse', 'hypre' 2199 2200 Level: intermediate 2201 2202 .seealso: PCMGGalerkinSetMatProductAlgorithm() 2203 2204 @*/ 2205 PetscErrorCode PCMGGalerkinGetMatProductAlgorithm(PC pc,const char *name[]) 2206 { 2207 PetscFunctionBegin; 2208 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 2209 CHKERRQ(PetscTryMethod(pc,"PCMGGalerkinGetMatProductAlgorithm_C",(PC,const char*[]),(pc,name))); 2210 PetscFunctionReturn(0); 2211 } 2212 2213 /*MC 2214 PCHYPRE - Allows you to use the matrix element based preconditioners in the LLNL package hypre 2215 2216 Options Database Keys: 2217 + -pc_hypre_type - One of euclid, pilut, parasails, boomeramg, ams, ads 2218 . -pc_hypre_boomeramg_nodal_coarsen <n> - where n is from 1 to 6 (see HYPRE_BOOMERAMGSetNodal()) 2219 . -pc_hypre_boomeramg_vec_interp_variant <v> - where v is from 1 to 3 (see HYPRE_BoomerAMGSetInterpVecVariant()) 2220 - Many others, run with -pc_type hypre -pc_hypre_type XXX -help to see options for the XXX preconditioner 2221 2222 Level: intermediate 2223 2224 Notes: 2225 Apart from pc_hypre_type (for which there is PCHYPRESetType()), 2226 the many hypre options can ONLY be set via the options database (e.g. the command line 2227 or with PetscOptionsSetValue(), there are no functions to set them) 2228 2229 The options -pc_hypre_boomeramg_max_iter and -pc_hypre_boomeramg_tol refer to the number of iterations 2230 (V-cycles) and tolerance that boomeramg does EACH time it is called. So for example, if 2231 -pc_hypre_boomeramg_max_iter is set to 2 then 2-V-cycles are being used to define the preconditioner 2232 (-pc_hypre_boomeramg_tol should be set to 0.0 - the default - to strictly use a fixed number of 2233 iterations per hypre call). -ksp_max_it and -ksp_rtol STILL determine the total number of iterations 2234 and tolerance for the Krylov solver. For example, if -pc_hypre_boomeramg_max_iter is 2 and -ksp_max_it is 10 2235 then AT MOST twenty V-cycles of boomeramg will be called. 2236 2237 Note that the option -pc_hypre_boomeramg_relax_type_all defaults to symmetric relaxation 2238 (symmetric-SOR/Jacobi), which is required for Krylov solvers like CG that expect symmetry. 2239 Otherwise, you may want to use -pc_hypre_boomeramg_relax_type_all SOR/Jacobi. 2240 If you wish to use BoomerAMG WITHOUT a Krylov method use -ksp_type richardson NOT -ksp_type preonly 2241 and use -ksp_max_it to control the number of V-cycles. 2242 (see the PETSc FAQ.html at the PETSc website under the Documentation tab). 2243 2244 2007-02-03 Using HYPRE-1.11.1b, the routine HYPRE_BoomerAMGSolveT and the option 2245 -pc_hypre_parasails_reuse were failing with SIGSEGV. Dalcin L. 2246 2247 MatSetNearNullSpace() - if you provide a near null space to your matrix it is ignored by hypre UNLESS you also use 2248 the following two options: 2249 2250 See PCPFMG for access to the hypre Struct PFMG solver 2251 2252 GPU Notes: 2253 To configure hypre BoomerAMG so that it can utilize NVIDIA GPUs run ./configure --download-hypre --with-cuda 2254 Then pass VECCUDA vectors and MATAIJCUSPARSE matrices to the solvers and PETSc will automatically utilize hypre's GPU solvers. 2255 2256 To configure hypre BoomerAMG so that it can utilize AMD GPUs run ./configure --download-hypre --with-hip 2257 Then pass VECHIP vectors to the solvers and PETSc will automatically utilize hypre's GPU solvers. 2258 2259 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 2260 PCHYPRESetType(), PCPFMG 2261 2262 M*/ 2263 2264 PETSC_EXTERN PetscErrorCode PCCreate_HYPRE(PC pc) 2265 { 2266 PC_HYPRE *jac; 2267 2268 PetscFunctionBegin; 2269 CHKERRQ(PetscNewLog(pc,&jac)); 2270 2271 pc->data = jac; 2272 pc->ops->reset = PCReset_HYPRE; 2273 pc->ops->destroy = PCDestroy_HYPRE; 2274 pc->ops->setfromoptions = PCSetFromOptions_HYPRE; 2275 pc->ops->setup = PCSetUp_HYPRE; 2276 pc->ops->apply = PCApply_HYPRE; 2277 jac->comm_hypre = MPI_COMM_NULL; 2278 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetType_C",PCHYPRESetType_HYPRE)); 2279 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPREGetType_C",PCHYPREGetType_HYPRE)); 2280 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCSetCoordinates_C",PCSetCoordinates_HYPRE)); 2281 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteGradient_C",PCHYPRESetDiscreteGradient_HYPRE)); 2282 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetDiscreteCurl_C",PCHYPRESetDiscreteCurl_HYPRE)); 2283 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetInterpolations_C",PCHYPRESetInterpolations_HYPRE)); 2284 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetEdgeConstantVectors_C",PCHYPRESetEdgeConstantVectors_HYPRE)); 2285 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCHYPRESetPoissonMatrix_C",PCHYPRESetPoissonMatrix_HYPRE)); 2286 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCMGGalerkinSetMatProductAlgorithm_C",PCMGGalerkinSetMatProductAlgorithm_HYPRE_BoomerAMG)); 2287 CHKERRQ(PetscObjectComposeFunction((PetscObject)pc,"PCMGGalerkinGetMatProductAlgorithm_C",PCMGGalerkinGetMatProductAlgorithm_HYPRE_BoomerAMG)); 2288 #if defined(PETSC_HAVE_HYPRE_DEVICE) 2289 #if defined(HYPRE_USING_HIP) 2290 CHKERRQ(PetscDeviceInitialize(PETSC_DEVICE_HIP)); 2291 #endif 2292 #if defined(HYPRE_USING_CUDA) 2293 CHKERRQ(PetscDeviceInitialize(PETSC_DEVICE_CUDA)); 2294 #endif 2295 #endif 2296 PetscFunctionReturn(0); 2297 } 2298 2299 /* ---------------------------------------------------------------------------------------------------------------------------------*/ 2300 2301 typedef struct { 2302 MPI_Comm hcomm; /* does not share comm with HYPRE_StructMatrix because need to create solver before getting matrix */ 2303 HYPRE_StructSolver hsolver; 2304 2305 /* keep copy of PFMG options used so may view them */ 2306 PetscInt its; 2307 double tol; 2308 PetscInt relax_type; 2309 PetscInt rap_type; 2310 PetscInt num_pre_relax,num_post_relax; 2311 PetscInt max_levels; 2312 } PC_PFMG; 2313 2314 PetscErrorCode PCDestroy_PFMG(PC pc) 2315 { 2316 PC_PFMG *ex = (PC_PFMG*) pc->data; 2317 2318 PetscFunctionBegin; 2319 if (ex->hsolver) PetscStackCallStandard(HYPRE_StructPFMGDestroy,ex->hsolver); 2320 CHKERRQ(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc),&ex->hcomm)); 2321 CHKERRQ(PetscFree(pc->data)); 2322 PetscFunctionReturn(0); 2323 } 2324 2325 static const char *PFMGRelaxType[] = {"Jacobi","Weighted-Jacobi","symmetric-Red/Black-Gauss-Seidel","Red/Black-Gauss-Seidel"}; 2326 static const char *PFMGRAPType[] = {"Galerkin","non-Galerkin"}; 2327 2328 PetscErrorCode PCView_PFMG(PC pc,PetscViewer viewer) 2329 { 2330 PetscBool iascii; 2331 PC_PFMG *ex = (PC_PFMG*) pc->data; 2332 2333 PetscFunctionBegin; 2334 CHKERRQ(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii)); 2335 if (iascii) { 2336 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE PFMG preconditioning\n")); 2337 CHKERRQ(PetscViewerASCIIPrintf(viewer," max iterations %d\n",ex->its)); 2338 CHKERRQ(PetscViewerASCIIPrintf(viewer," tolerance %g\n",ex->tol)); 2339 CHKERRQ(PetscViewerASCIIPrintf(viewer," relax type %s\n",PFMGRelaxType[ex->relax_type])); 2340 CHKERRQ(PetscViewerASCIIPrintf(viewer," RAP type %s\n",PFMGRAPType[ex->rap_type])); 2341 CHKERRQ(PetscViewerASCIIPrintf(viewer," number pre-relax %d post-relax %d\n",ex->num_pre_relax,ex->num_post_relax)); 2342 CHKERRQ(PetscViewerASCIIPrintf(viewer," max levels %d\n",ex->max_levels)); 2343 } 2344 PetscFunctionReturn(0); 2345 } 2346 2347 PetscErrorCode PCSetFromOptions_PFMG(PetscOptionItems *PetscOptionsObject,PC pc) 2348 { 2349 PC_PFMG *ex = (PC_PFMG*) pc->data; 2350 PetscBool flg = PETSC_FALSE; 2351 2352 PetscFunctionBegin; 2353 CHKERRQ(PetscOptionsHead(PetscOptionsObject,"PFMG options")); 2354 CHKERRQ(PetscOptionsBool("-pc_pfmg_print_statistics","Print statistics","HYPRE_StructPFMGSetPrintLevel",flg,&flg,NULL)); 2355 if (flg) { 2356 PetscStackCallStandard(HYPRE_StructPFMGSetPrintLevel,ex->hsolver,3); 2357 } 2358 CHKERRQ(PetscOptionsInt("-pc_pfmg_its","Number of iterations of PFMG to use as preconditioner","HYPRE_StructPFMGSetMaxIter",ex->its,&ex->its,NULL)); 2359 PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,ex->hsolver,ex->its); 2360 CHKERRQ(PetscOptionsInt("-pc_pfmg_num_pre_relax","Number of smoothing steps before coarse grid","HYPRE_StructPFMGSetNumPreRelax",ex->num_pre_relax,&ex->num_pre_relax,NULL)); 2361 PetscStackCallStandard(HYPRE_StructPFMGSetNumPreRelax,ex->hsolver,ex->num_pre_relax); 2362 CHKERRQ(PetscOptionsInt("-pc_pfmg_num_post_relax","Number of smoothing steps after coarse grid","HYPRE_StructPFMGSetNumPostRelax",ex->num_post_relax,&ex->num_post_relax,NULL)); 2363 PetscStackCallStandard(HYPRE_StructPFMGSetNumPostRelax,ex->hsolver,ex->num_post_relax); 2364 2365 CHKERRQ(PetscOptionsInt("-pc_pfmg_max_levels","Max Levels for MG hierarchy","HYPRE_StructPFMGSetMaxLevels",ex->max_levels,&ex->max_levels,NULL)); 2366 PetscStackCallStandard(HYPRE_StructPFMGSetMaxLevels,ex->hsolver,ex->max_levels); 2367 2368 CHKERRQ(PetscOptionsReal("-pc_pfmg_tol","Tolerance of PFMG","HYPRE_StructPFMGSetTol",ex->tol,&ex->tol,NULL)); 2369 PetscStackCallStandard(HYPRE_StructPFMGSetTol,ex->hsolver,ex->tol); 2370 CHKERRQ(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)); 2371 PetscStackCallStandard(HYPRE_StructPFMGSetRelaxType,ex->hsolver, ex->relax_type); 2372 CHKERRQ(PetscOptionsEList("-pc_pfmg_rap_type","RAP type","HYPRE_StructPFMGSetRAPType",PFMGRAPType,ALEN(PFMGRAPType),PFMGRAPType[ex->rap_type],&ex->rap_type,NULL)); 2373 PetscStackCallStandard(HYPRE_StructPFMGSetRAPType,ex->hsolver, ex->rap_type); 2374 CHKERRQ(PetscOptionsTail()); 2375 PetscFunctionReturn(0); 2376 } 2377 2378 PetscErrorCode PCApply_PFMG(PC pc,Vec x,Vec y) 2379 { 2380 PC_PFMG *ex = (PC_PFMG*) pc->data; 2381 PetscScalar *yy; 2382 const PetscScalar *xx; 2383 PetscInt ilower[3],iupper[3]; 2384 HYPRE_Int hlower[3],hupper[3]; 2385 Mat_HYPREStruct *mx = (Mat_HYPREStruct*)(pc->pmat->data); 2386 2387 PetscFunctionBegin; 2388 CHKERRQ(PetscCitationsRegister(hypreCitation,&cite)); 2389 CHKERRQ(DMDAGetCorners(mx->da,&ilower[0],&ilower[1],&ilower[2],&iupper[0],&iupper[1],&iupper[2])); 2390 /* when HYPRE_MIXEDINT is defined, sizeof(HYPRE_Int) == 32 */ 2391 iupper[0] += ilower[0] - 1; 2392 iupper[1] += ilower[1] - 1; 2393 iupper[2] += ilower[2] - 1; 2394 hlower[0] = (HYPRE_Int)ilower[0]; 2395 hlower[1] = (HYPRE_Int)ilower[1]; 2396 hlower[2] = (HYPRE_Int)ilower[2]; 2397 hupper[0] = (HYPRE_Int)iupper[0]; 2398 hupper[1] = (HYPRE_Int)iupper[1]; 2399 hupper[2] = (HYPRE_Int)iupper[2]; 2400 2401 /* copy x values over to hypre */ 2402 PetscStackCallStandard(HYPRE_StructVectorSetConstantValues,mx->hb,0.0); 2403 CHKERRQ(VecGetArrayRead(x,&xx)); 2404 PetscStackCallStandard(HYPRE_StructVectorSetBoxValues,mx->hb,hlower,hupper,(HYPRE_Complex*)xx); 2405 CHKERRQ(VecRestoreArrayRead(x,&xx)); 2406 PetscStackCallStandard(HYPRE_StructVectorAssemble,mx->hb); 2407 PetscStackCallStandard(HYPRE_StructPFMGSolve,ex->hsolver,mx->hmat,mx->hb,mx->hx); 2408 2409 /* copy solution values back to PETSc */ 2410 CHKERRQ(VecGetArray(y,&yy)); 2411 PetscStackCallStandard(HYPRE_StructVectorGetBoxValues,mx->hx,hlower,hupper,(HYPRE_Complex*)yy); 2412 CHKERRQ(VecRestoreArray(y,&yy)); 2413 PetscFunctionReturn(0); 2414 } 2415 2416 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) 2417 { 2418 PC_PFMG *jac = (PC_PFMG*)pc->data; 2419 HYPRE_Int oits; 2420 2421 PetscFunctionBegin; 2422 CHKERRQ(PetscCitationsRegister(hypreCitation,&cite)); 2423 PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,jac->hsolver,its*jac->its); 2424 PetscStackCallStandard(HYPRE_StructPFMGSetTol,jac->hsolver,rtol); 2425 2426 CHKERRQ(PCApply_PFMG(pc,b,y)); 2427 PetscStackCallStandard(HYPRE_StructPFMGGetNumIterations,jac->hsolver,&oits); 2428 *outits = oits; 2429 if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS; 2430 else *reason = PCRICHARDSON_CONVERGED_RTOL; 2431 PetscStackCallStandard(HYPRE_StructPFMGSetTol,jac->hsolver,jac->tol); 2432 PetscStackCallStandard(HYPRE_StructPFMGSetMaxIter,jac->hsolver,jac->its); 2433 PetscFunctionReturn(0); 2434 } 2435 2436 PetscErrorCode PCSetUp_PFMG(PC pc) 2437 { 2438 PC_PFMG *ex = (PC_PFMG*) pc->data; 2439 Mat_HYPREStruct *mx = (Mat_HYPREStruct*)(pc->pmat->data); 2440 PetscBool flg; 2441 2442 PetscFunctionBegin; 2443 CHKERRQ(PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRESTRUCT,&flg)); 2444 PetscCheck(flg,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Must use MATHYPRESTRUCT with this preconditioner"); 2445 2446 /* create the hypre solver object and set its information */ 2447 if (ex->hsolver) PetscStackCallStandard(HYPRE_StructPFMGDestroy,ex->hsolver); 2448 PetscStackCallStandard(HYPRE_StructPFMGCreate,ex->hcomm,&ex->hsolver); 2449 PetscStackCallStandard(HYPRE_StructPFMGSetup,ex->hsolver,mx->hmat,mx->hb,mx->hx); 2450 PetscStackCallStandard(HYPRE_StructPFMGSetZeroGuess,ex->hsolver); 2451 PetscFunctionReturn(0); 2452 } 2453 2454 /*MC 2455 PCPFMG - the hypre PFMG multigrid solver 2456 2457 Level: advanced 2458 2459 Options Database: 2460 + -pc_pfmg_its <its> - number of iterations of PFMG to use as preconditioner 2461 . -pc_pfmg_num_pre_relax <steps> - number of smoothing steps before coarse grid solve 2462 . -pc_pfmg_num_post_relax <steps> - number of smoothing steps after coarse grid solve 2463 . -pc_pfmg_tol <tol> - tolerance of PFMG 2464 . -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 2465 - -pc_pfmg_rap_type - type of coarse matrix generation, one of Galerkin,non-Galerkin 2466 2467 Notes: 2468 This is for CELL-centered descretizations 2469 2470 This must be used with the MATHYPRESTRUCT matrix type. 2471 This is less general than in hypre, it supports only one block per process defined by a PETSc DMDA. 2472 2473 .seealso: PCMG, MATHYPRESTRUCT 2474 M*/ 2475 2476 PETSC_EXTERN PetscErrorCode PCCreate_PFMG(PC pc) 2477 { 2478 PC_PFMG *ex; 2479 2480 PetscFunctionBegin; 2481 CHKERRQ(PetscNew(&ex)); \ 2482 pc->data = ex; 2483 2484 ex->its = 1; 2485 ex->tol = 1.e-8; 2486 ex->relax_type = 1; 2487 ex->rap_type = 0; 2488 ex->num_pre_relax = 1; 2489 ex->num_post_relax = 1; 2490 ex->max_levels = 0; 2491 2492 pc->ops->setfromoptions = PCSetFromOptions_PFMG; 2493 pc->ops->view = PCView_PFMG; 2494 pc->ops->destroy = PCDestroy_PFMG; 2495 pc->ops->apply = PCApply_PFMG; 2496 pc->ops->applyrichardson = PCApplyRichardson_PFMG; 2497 pc->ops->setup = PCSetUp_PFMG; 2498 2499 CHKERRQ(PetscCommGetComm(PetscObjectComm((PetscObject)pc),&ex->hcomm)); 2500 PetscStackCallStandard(HYPRE_StructPFMGCreate,ex->hcomm,&ex->hsolver); 2501 PetscFunctionReturn(0); 2502 } 2503 2504 /* ---------------------------------------------------------------------------------------------------------------------------------------------------*/ 2505 2506 /* we know we are working with a HYPRE_SStructMatrix */ 2507 typedef struct { 2508 MPI_Comm hcomm; /* does not share comm with HYPRE_SStructMatrix because need to create solver before getting matrix */ 2509 HYPRE_SStructSolver ss_solver; 2510 2511 /* keep copy of SYSPFMG options used so may view them */ 2512 PetscInt its; 2513 double tol; 2514 PetscInt relax_type; 2515 PetscInt num_pre_relax,num_post_relax; 2516 } PC_SysPFMG; 2517 2518 PetscErrorCode PCDestroy_SysPFMG(PC pc) 2519 { 2520 PC_SysPFMG *ex = (PC_SysPFMG*) pc->data; 2521 2522 PetscFunctionBegin; 2523 if (ex->ss_solver) PetscStackCallStandard(HYPRE_SStructSysPFMGDestroy,ex->ss_solver); 2524 CHKERRQ(PetscCommRestoreComm(PetscObjectComm((PetscObject)pc),&ex->hcomm)); 2525 CHKERRQ(PetscFree(pc->data)); 2526 PetscFunctionReturn(0); 2527 } 2528 2529 static const char *SysPFMGRelaxType[] = {"Weighted-Jacobi","Red/Black-Gauss-Seidel"}; 2530 2531 PetscErrorCode PCView_SysPFMG(PC pc,PetscViewer viewer) 2532 { 2533 PetscBool iascii; 2534 PC_SysPFMG *ex = (PC_SysPFMG*) pc->data; 2535 2536 PetscFunctionBegin; 2537 CHKERRQ(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii)); 2538 if (iascii) { 2539 CHKERRQ(PetscViewerASCIIPrintf(viewer," HYPRE SysPFMG preconditioning\n")); 2540 CHKERRQ(PetscViewerASCIIPrintf(viewer," max iterations %d\n",ex->its)); 2541 CHKERRQ(PetscViewerASCIIPrintf(viewer," tolerance %g\n",ex->tol)); 2542 CHKERRQ(PetscViewerASCIIPrintf(viewer," relax type %s\n",PFMGRelaxType[ex->relax_type])); 2543 CHKERRQ(PetscViewerASCIIPrintf(viewer," number pre-relax %d post-relax %d\n",ex->num_pre_relax,ex->num_post_relax)); 2544 } 2545 PetscFunctionReturn(0); 2546 } 2547 2548 PetscErrorCode PCSetFromOptions_SysPFMG(PetscOptionItems *PetscOptionsObject,PC pc) 2549 { 2550 PC_SysPFMG *ex = (PC_SysPFMG*) pc->data; 2551 PetscBool flg = PETSC_FALSE; 2552 2553 PetscFunctionBegin; 2554 CHKERRQ(PetscOptionsHead(PetscOptionsObject,"SysPFMG options")); 2555 CHKERRQ(PetscOptionsBool("-pc_syspfmg_print_statistics","Print statistics","HYPRE_SStructSysPFMGSetPrintLevel",flg,&flg,NULL)); 2556 if (flg) { 2557 PetscStackCallStandard(HYPRE_SStructSysPFMGSetPrintLevel,ex->ss_solver,3); 2558 } 2559 CHKERRQ(PetscOptionsInt("-pc_syspfmg_its","Number of iterations of SysPFMG to use as preconditioner","HYPRE_SStructSysPFMGSetMaxIter",ex->its,&ex->its,NULL)); 2560 PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,ex->ss_solver,ex->its); 2561 CHKERRQ(PetscOptionsInt("-pc_syspfmg_num_pre_relax","Number of smoothing steps before coarse grid","HYPRE_SStructSysPFMGSetNumPreRelax",ex->num_pre_relax,&ex->num_pre_relax,NULL)); 2562 PetscStackCallStandard(HYPRE_SStructSysPFMGSetNumPreRelax,ex->ss_solver,ex->num_pre_relax); 2563 CHKERRQ(PetscOptionsInt("-pc_syspfmg_num_post_relax","Number of smoothing steps after coarse grid","HYPRE_SStructSysPFMGSetNumPostRelax",ex->num_post_relax,&ex->num_post_relax,NULL)); 2564 PetscStackCallStandard(HYPRE_SStructSysPFMGSetNumPostRelax,ex->ss_solver,ex->num_post_relax); 2565 2566 CHKERRQ(PetscOptionsReal("-pc_syspfmg_tol","Tolerance of SysPFMG","HYPRE_SStructSysPFMGSetTol",ex->tol,&ex->tol,NULL)); 2567 PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,ex->ss_solver,ex->tol); 2568 CHKERRQ(PetscOptionsEList("-pc_syspfmg_relax_type","Relax type for the up and down cycles","HYPRE_SStructSysPFMGSetRelaxType",SysPFMGRelaxType,ALEN(SysPFMGRelaxType),SysPFMGRelaxType[ex->relax_type],&ex->relax_type,NULL)); 2569 PetscStackCallStandard(HYPRE_SStructSysPFMGSetRelaxType,ex->ss_solver, ex->relax_type); 2570 CHKERRQ(PetscOptionsTail()); 2571 PetscFunctionReturn(0); 2572 } 2573 2574 PetscErrorCode PCApply_SysPFMG(PC pc,Vec x,Vec y) 2575 { 2576 PC_SysPFMG *ex = (PC_SysPFMG*) pc->data; 2577 PetscScalar *yy; 2578 const PetscScalar *xx; 2579 PetscInt ilower[3],iupper[3]; 2580 HYPRE_Int hlower[3],hupper[3]; 2581 Mat_HYPRESStruct *mx = (Mat_HYPRESStruct*)(pc->pmat->data); 2582 PetscInt ordering= mx->dofs_order; 2583 PetscInt nvars = mx->nvars; 2584 PetscInt part = 0; 2585 PetscInt size; 2586 PetscInt i; 2587 2588 PetscFunctionBegin; 2589 CHKERRQ(PetscCitationsRegister(hypreCitation,&cite)); 2590 CHKERRQ(DMDAGetCorners(mx->da,&ilower[0],&ilower[1],&ilower[2],&iupper[0],&iupper[1],&iupper[2])); 2591 /* when HYPRE_MIXEDINT is defined, sizeof(HYPRE_Int) == 32 */ 2592 iupper[0] += ilower[0] - 1; 2593 iupper[1] += ilower[1] - 1; 2594 iupper[2] += ilower[2] - 1; 2595 hlower[0] = (HYPRE_Int)ilower[0]; 2596 hlower[1] = (HYPRE_Int)ilower[1]; 2597 hlower[2] = (HYPRE_Int)ilower[2]; 2598 hupper[0] = (HYPRE_Int)iupper[0]; 2599 hupper[1] = (HYPRE_Int)iupper[1]; 2600 hupper[2] = (HYPRE_Int)iupper[2]; 2601 2602 size = 1; 2603 for (i= 0; i< 3; i++) size *= (iupper[i]-ilower[i]+1); 2604 2605 /* copy x values over to hypre for variable ordering */ 2606 if (ordering) { 2607 PetscStackCallStandard(HYPRE_SStructVectorSetConstantValues,mx->ss_b,0.0); 2608 CHKERRQ(VecGetArrayRead(x,&xx)); 2609 for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorSetBoxValues,mx->ss_b,part,hlower,hupper,i,(HYPRE_Complex*)(xx+(size*i))); 2610 CHKERRQ(VecRestoreArrayRead(x,&xx)); 2611 PetscStackCallStandard(HYPRE_SStructVectorAssemble,mx->ss_b); 2612 PetscStackCallStandard(HYPRE_SStructMatrixMatvec,1.0,mx->ss_mat,mx->ss_b,0.0,mx->ss_x); 2613 PetscStackCallStandard(HYPRE_SStructSysPFMGSolve,ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x); 2614 2615 /* copy solution values back to PETSc */ 2616 CHKERRQ(VecGetArray(y,&yy)); 2617 for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorGetBoxValues,mx->ss_x,part,hlower,hupper,i,(HYPRE_Complex*)(yy+(size*i))); 2618 CHKERRQ(VecRestoreArray(y,&yy)); 2619 } else { /* nodal ordering must be mapped to variable ordering for sys_pfmg */ 2620 PetscScalar *z; 2621 PetscInt j, k; 2622 2623 CHKERRQ(PetscMalloc1(nvars*size,&z)); 2624 PetscStackCallStandard(HYPRE_SStructVectorSetConstantValues,mx->ss_b,0.0); 2625 CHKERRQ(VecGetArrayRead(x,&xx)); 2626 2627 /* transform nodal to hypre's variable ordering for sys_pfmg */ 2628 for (i= 0; i< size; i++) { 2629 k= i*nvars; 2630 for (j= 0; j< nvars; j++) z[j*size+i]= xx[k+j]; 2631 } 2632 for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorSetBoxValues,mx->ss_b,part,hlower,hupper,i,(HYPRE_Complex*)(z+(size*i))); 2633 CHKERRQ(VecRestoreArrayRead(x,&xx)); 2634 PetscStackCallStandard(HYPRE_SStructVectorAssemble,mx->ss_b); 2635 PetscStackCallStandard(HYPRE_SStructSysPFMGSolve,ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x); 2636 2637 /* copy solution values back to PETSc */ 2638 CHKERRQ(VecGetArray(y,&yy)); 2639 for (i= 0; i< nvars; i++) PetscStackCallStandard(HYPRE_SStructVectorGetBoxValues,mx->ss_x,part,hlower,hupper,i,(HYPRE_Complex*)(z+(size*i))); 2640 /* transform hypre's variable ordering for sys_pfmg to nodal ordering */ 2641 for (i= 0; i< size; i++) { 2642 k= i*nvars; 2643 for (j= 0; j< nvars; j++) yy[k+j]= z[j*size+i]; 2644 } 2645 CHKERRQ(VecRestoreArray(y,&yy)); 2646 CHKERRQ(PetscFree(z)); 2647 } 2648 PetscFunctionReturn(0); 2649 } 2650 2651 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) 2652 { 2653 PC_SysPFMG *jac = (PC_SysPFMG*)pc->data; 2654 HYPRE_Int oits; 2655 2656 PetscFunctionBegin; 2657 CHKERRQ(PetscCitationsRegister(hypreCitation,&cite)); 2658 PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,jac->ss_solver,its*jac->its); 2659 PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,jac->ss_solver,rtol); 2660 CHKERRQ(PCApply_SysPFMG(pc,b,y)); 2661 PetscStackCallStandard(HYPRE_SStructSysPFMGGetNumIterations,jac->ss_solver,&oits); 2662 *outits = oits; 2663 if (oits == its) *reason = PCRICHARDSON_CONVERGED_ITS; 2664 else *reason = PCRICHARDSON_CONVERGED_RTOL; 2665 PetscStackCallStandard(HYPRE_SStructSysPFMGSetTol,jac->ss_solver,jac->tol); 2666 PetscStackCallStandard(HYPRE_SStructSysPFMGSetMaxIter,jac->ss_solver,jac->its); 2667 PetscFunctionReturn(0); 2668 } 2669 2670 PetscErrorCode PCSetUp_SysPFMG(PC pc) 2671 { 2672 PC_SysPFMG *ex = (PC_SysPFMG*) pc->data; 2673 Mat_HYPRESStruct *mx = (Mat_HYPRESStruct*)(pc->pmat->data); 2674 PetscBool flg; 2675 2676 PetscFunctionBegin; 2677 CHKERRQ(PetscObjectTypeCompare((PetscObject)pc->pmat,MATHYPRESSTRUCT,&flg)); 2678 PetscCheck(flg,PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_INCOMP,"Must use MATHYPRESSTRUCT with this preconditioner"); 2679 2680 /* create the hypre sstruct solver object and set its information */ 2681 if (ex->ss_solver) PetscStackCallStandard(HYPRE_SStructSysPFMGDestroy,ex->ss_solver); 2682 PetscStackCallStandard(HYPRE_SStructSysPFMGCreate,ex->hcomm,&ex->ss_solver); 2683 PetscStackCallStandard(HYPRE_SStructSysPFMGSetZeroGuess,ex->ss_solver); 2684 PetscStackCallStandard(HYPRE_SStructSysPFMGSetup,ex->ss_solver,mx->ss_mat,mx->ss_b,mx->ss_x); 2685 PetscFunctionReturn(0); 2686 } 2687 2688 /*MC 2689 PCSysPFMG - the hypre SysPFMG multigrid solver 2690 2691 Level: advanced 2692 2693 Options Database: 2694 + -pc_syspfmg_its <its> - number of iterations of SysPFMG to use as preconditioner 2695 . -pc_syspfmg_num_pre_relax <steps> - number of smoothing steps before coarse grid 2696 . -pc_syspfmg_num_post_relax <steps> - number of smoothing steps after coarse grid 2697 . -pc_syspfmg_tol <tol> - tolerance of SysPFMG 2698 - -pc_syspfmg_relax_type <Weighted-Jacobi,Red/Black-Gauss-Seidel> - relaxation type for the up and down cycles 2699 2700 Notes: 2701 This is for CELL-centered descretizations 2702 2703 This must be used with the MATHYPRESSTRUCT matrix type. 2704 This is less general than in hypre, it supports only one part, and one block per process defined by a PETSc DMDA. 2705 Also, only cell-centered variables. 2706 2707 .seealso: PCMG, MATHYPRESSTRUCT 2708 M*/ 2709 2710 PETSC_EXTERN PetscErrorCode PCCreate_SysPFMG(PC pc) 2711 { 2712 PC_SysPFMG *ex; 2713 2714 PetscFunctionBegin; 2715 CHKERRQ(PetscNew(&ex)); \ 2716 pc->data = ex; 2717 2718 ex->its = 1; 2719 ex->tol = 1.e-8; 2720 ex->relax_type = 1; 2721 ex->num_pre_relax = 1; 2722 ex->num_post_relax = 1; 2723 2724 pc->ops->setfromoptions = PCSetFromOptions_SysPFMG; 2725 pc->ops->view = PCView_SysPFMG; 2726 pc->ops->destroy = PCDestroy_SysPFMG; 2727 pc->ops->apply = PCApply_SysPFMG; 2728 pc->ops->applyrichardson = PCApplyRichardson_SysPFMG; 2729 pc->ops->setup = PCSetUp_SysPFMG; 2730 2731 CHKERRQ(PetscCommGetComm(PetscObjectComm((PetscObject)pc),&ex->hcomm)); 2732 PetscStackCallStandard(HYPRE_SStructSysPFMGCreate,ex->hcomm,&ex->ss_solver); 2733 PetscFunctionReturn(0); 2734 } 2735