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