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