1 /* 2 GAMG geometric-algebric multiogrid PC - Mark Adams 2011 3 */ 4 #include <private/pcimpl.h> /*I "petscpc.h" I*/ 5 #include <../src/ksp/pc/impls/mg/mgimpl.h> /*I "petscpcmg.h" I*/ 6 #include <../src/mat/impls/aij/seq/aij.h> 7 #include <../src/mat/impls/aij/mpi/mpiaij.h> 8 #include <assert.h> 9 10 /* Private context for the GAMG preconditioner */ 11 typedef struct gamg_TAG { 12 PetscInt m_dim; 13 PetscInt m_Nlevels; 14 PetscInt m_data_sz; 15 PetscReal *m_data; /* blocked vector of vertex data on fine grid (coordinates) */ 16 } PC_GAMG; 17 18 /* -----------------------------------------------------------------------------*/ 19 #undef __FUNCT__ 20 #define __FUNCT__ "PCReset_GAMG" 21 PetscErrorCode PCReset_GAMG(PC pc) 22 { 23 PetscErrorCode ierr; 24 PC_MG *mg = (PC_MG*)pc->data; 25 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 26 27 PetscFunctionBegin; 28 ierr = PetscFree(pc_gamg->m_data);CHKERRQ(ierr); 29 PetscFunctionReturn(0); 30 } 31 32 /* -------------------------------------------------------------------------- */ 33 /* 34 PCSetCoordinates_GAMG 35 36 Input Parameter: 37 . pc - the preconditioner context 38 */ 39 EXTERN_C_BEGIN 40 #undef __FUNCT__ 41 #define __FUNCT__ "PCSetCoordinates_GAMG" 42 PetscErrorCode PCSetCoordinates_GAMG(PC pc, const int ndm,PetscReal *coords ) 43 { 44 PC_MG *mg = (PC_MG*)pc->data; 45 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 46 PetscErrorCode ierr; 47 PetscInt bs, my0, tt; 48 Mat mat = pc->pmat; 49 PetscInt arrsz; 50 51 PetscFunctionBegin; 52 ierr = MatGetBlockSize( mat, &bs ); CHKERRQ( ierr ); 53 ierr = MatGetOwnershipRange( mat, &my0, &tt ); CHKERRQ(ierr); 54 arrsz = (tt-my0)/bs*ndm; 55 56 // put coordinates 57 if (!pc_gamg->m_data || (pc_gamg->m_data_sz != arrsz)) { 58 ierr = PetscFree( pc_gamg->m_data ); CHKERRQ(ierr); 59 ierr = PetscMalloc(arrsz*sizeof(double), &pc_gamg->m_data ); CHKERRQ(ierr); 60 } 61 62 /* copy data in */ 63 for(tt=0;tt<arrsz;tt++){ 64 pc_gamg->m_data[tt] = coords[tt]; 65 } 66 pc_gamg->m_data_sz = arrsz; 67 pc_gamg->m_dim = ndm; 68 69 PetscFunctionReturn(0); 70 } 71 EXTERN_C_END 72 73 /* -------------------------------------------------------------------------- */ 74 /* 75 partitionLevel 76 77 Input Parameter: 78 . a_Amat_fine - matrix on this fine (k) level 79 . a_dime - 2 or 3 80 In/Output Parameter: 81 . a_P_inout - prolongation operator to the next level (k-1) 82 . a_coarse_crds - coordinates that need to be moved 83 . a_active_proc - number of active procs 84 Output Parameter: 85 . a_Amat_crs - coarse matrix that is created (k-1) 86 */ 87 #undef __FUNCT__ 88 #define __FUNCT__ "partitionLevel" 89 PetscErrorCode partitionLevel( Mat a_Amat_fine, 90 PetscInt a_dim, 91 Mat *a_P_inout, 92 PetscReal **a_coarse_crds, 93 PetscMPIInt *a_active_proc, 94 Mat *a_Amat_crs 95 ) 96 { 97 PetscErrorCode ierr; 98 Mat Amat, Pnew, Pold = *a_P_inout; 99 IS new_indices,isnum; 100 MPI_Comm wcomm = ((PetscObject)a_Amat_fine)->comm; 101 PetscMPIInt nactive,mype,npe; 102 PetscInt Istart,Iend,Istart0,Iend0,ncrs0,ncrs_new,bs=1; /* bs ??? */ 103 104 PetscFunctionBegin; 105 ierr = MPI_Comm_rank(wcomm,&mype);CHKERRQ(ierr); 106 ierr = MPI_Comm_size(wcomm,&npe);CHKERRQ(ierr); 107 /* RAP */ 108 ierr = MatPtAP( a_Amat_fine, Pold, MAT_INITIAL_MATRIX, 2.0, &Amat ); CHKERRQ(ierr); 109 ierr = MatGetOwnershipRange( Amat, &Istart0, &Iend0 ); CHKERRQ(ierr); /* x2 size of 'a_coarse_crds' */ 110 ncrs0 = (Iend0 - Istart0)/bs; 111 112 /* Repartition Amat_{k} and move colums of P^{k}_{k-1} and coordinates accordingly */ 113 { 114 PetscInt neq,N,counts[npe]; 115 IS isnewproc; 116 PetscMPIInt factstart,fact,new_npe,targ_npe; 117 118 ierr = MatGetSize( Amat, &neq, &N );CHKERRQ(ierr); 119 #define MIN_EQ_PROC 100 120 nactive = *a_active_proc; 121 targ_npe = neq/MIN_EQ_PROC; /* hardwire min. number of eq/proc */ 122 if( targ_npe == 0 || neq < 1000 ) new_npe = 1; /* chop coarsest grid */ 123 else if (targ_npe > nactive ) new_npe = nactive; /* no change */ 124 else { 125 new_npe = -9999; 126 factstart = nactive; 127 for(fact=factstart;fact>0;fact--){ /* try to find a better number of procs */ 128 if( nactive%fact==0 && neq/(nactive/fact) > MIN_EQ_PROC ) { 129 new_npe = nactive/fact; 130 } 131 } 132 assert(new_npe != -9999); 133 } 134 *a_active_proc = new_npe; /* output for next time */ 135 136 { /* partition: get 'isnewproc' */ 137 MatPartitioning mpart; 138 ierr = MatPartitioningCreate( wcomm, &mpart ); CHKERRQ(ierr); 139 ierr = MatPartitioningSetAdjacency( mpart, Amat ); CHKERRQ(ierr); 140 ierr = MatPartitioningSetFromOptions( mpart ); CHKERRQ(ierr); 141 ierr = MatPartitioningApply( mpart, &isnewproc ); CHKERRQ(ierr); 142 ierr = MatPartitioningDestroy( &mpart ); CHKERRQ(ierr); 143 } 144 /* 145 Create an index set from the isnewproc index set to indicate the mapping TO 146 */ 147 ierr = ISPartitioningToNumbering( isnewproc, &isnum ); CHKERRQ(ierr); 148 /* 149 Determine how many elements are assigned to each processor 150 */ 151 ierr = ISPartitioningCount( isnewproc, npe, counts ); CHKERRQ(ierr); 152 ierr = ISDestroy( &isnewproc ); CHKERRQ(ierr); 153 ncrs_new = counts[mype]; 154 } 155 { /* Create a vector to contain the newly ordered element information */ 156 const PetscInt *idx; 157 PetscInt i,j,k; 158 IS isscat; 159 PetscScalar *array; 160 Vec src_crd, dest_crd; 161 PetscReal *coords = *a_coarse_crds; 162 VecScatter vecscat; 163 164 ierr = VecCreate( wcomm, &dest_crd ); 165 ierr = VecSetSizes( dest_crd, a_dim*ncrs_new, PETSC_DECIDE ); CHKERRQ(ierr); 166 ierr = VecSetFromOptions( dest_crd ); CHKERRQ(ierr); /*funny vector-get global options?*/ 167 /* 168 There are three data items per cell (element), the integer vertex numbers of its three 169 coordinates (we convert to double to use the scatter) (one can think 170 of the vectors of having a block size of 3, then there is one index in idx[] for each element) 171 */ 172 { 173 PetscInt tidx[ncrs0*a_dim]; 174 ierr = ISGetIndices( isnum, &idx ); CHKERRQ(ierr); 175 for (i=0,j=0; i<ncrs0; i++) { 176 for (k=0; k<a_dim; k++) tidx[j++] = idx[i]*a_dim + k; 177 } 178 ierr = ISCreateGeneral( wcomm, a_dim*ncrs0, tidx, PETSC_COPY_VALUES, &isscat ); 179 CHKERRQ(ierr); 180 ierr = ISRestoreIndices( isnum, &idx ); CHKERRQ(ierr); 181 } 182 /* 183 Create a vector to contain the original vertex information for each element 184 */ 185 ierr = VecCreateSeq( PETSC_COMM_SELF, a_dim*ncrs0, &src_crd ); CHKERRQ(ierr); 186 ierr = VecGetArray( src_crd, &array ); CHKERRQ(ierr); 187 for (i=0; i<a_dim*ncrs0; i++) array[i] = coords[i]; 188 ierr = VecRestoreArray( src_crd, &array ); CHKERRQ(ierr); 189 /* 190 Scatter the element vertex information (still in the original vertex ordering) 191 to the correct processor 192 */ 193 ierr = VecScatterCreate( src_crd, PETSC_NULL, dest_crd, isscat, &vecscat); 194 CHKERRQ(ierr); 195 ierr = ISDestroy( &isscat ); CHKERRQ(ierr); 196 ierr = VecScatterBegin(vecscat,src_crd,dest_crd,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 197 ierr = VecScatterEnd(vecscat,src_crd,dest_crd,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 198 ierr = VecScatterDestroy( &vecscat ); CHKERRQ(ierr); 199 ierr = VecDestroy( &src_crd ); CHKERRQ(ierr); 200 /* 201 Put the element vertex data into a new allocation of the gdata->ele 202 */ 203 ierr = PetscFree( *a_coarse_crds ); CHKERRQ(ierr); 204 ierr = PetscMalloc( a_dim*ncrs_new*sizeof(PetscReal), a_coarse_crds ); CHKERRQ(ierr); 205 coords = *a_coarse_crds; /* convient */ 206 ierr = VecGetArray( dest_crd, &array ); CHKERRQ(ierr); 207 for (i=0; i<a_dim*ncrs_new; i++) coords[i] = PetscRealPart(array[i]); 208 ierr = VecRestoreArray( dest_crd, &array ); CHKERRQ(ierr); 209 ierr = VecDestroy( &dest_crd ); CHKERRQ(ierr); 210 } 211 /* 212 Invert for MatGetSubMatrix 213 */ 214 ierr = ISInvertPermutation( isnum, ncrs_new, &new_indices ); CHKERRQ(ierr); 215 ierr = ISSort( new_indices ); CHKERRQ(ierr); /* is this needed? */ 216 ierr = ISDestroy( &isnum ); CHKERRQ(ierr); 217 /* A_crs output */ 218 ierr = MatGetSubMatrix( Amat, new_indices, new_indices, MAT_INITIAL_MATRIX, a_Amat_crs ); 219 CHKERRQ(ierr); 220 ierr = MatDestroy( &Amat ); CHKERRQ(ierr); 221 Amat = *a_Amat_crs; 222 /* prolongator */ 223 ierr = MatGetOwnershipRange( Pold, &Istart, &Iend ); CHKERRQ(ierr); 224 { 225 IS findices; 226 ierr = ISCreateStride(wcomm,Iend-Istart,Istart,1,&findices); CHKERRQ(ierr); 227 ierr = MatGetSubMatrix( Pold, findices, new_indices, MAT_INITIAL_MATRIX, &Pnew ); 228 CHKERRQ(ierr); 229 ierr = ISDestroy( &findices ); CHKERRQ(ierr); 230 } 231 ierr = MatDestroy( a_P_inout ); CHKERRQ(ierr); 232 *a_P_inout = Pnew; /* output */ 233 234 ierr = ISDestroy( &new_indices ); CHKERRQ(ierr); 235 236 PetscFunctionReturn(0); 237 } 238 239 /* -------------------------------------------------------------------------- */ 240 /* 241 PCSetUp_GAMG - Prepares for the use of the GAMG preconditioner 242 by setting data structures and options. 243 244 Input Parameter: 245 . pc - the preconditioner context 246 247 Application Interface Routine: PCSetUp() 248 249 Notes: 250 The interface routine PCSetUp() is not usually called directly by 251 the user, but instead is called by PCApply() if necessary. 252 */ 253 extern PetscErrorCode PCSetFromOptions_MG(PC); 254 extern PetscErrorCode PCReset_MG(PC); 255 extern PetscErrorCode createProlongation( Mat, PetscReal [], const PetscInt, 256 Mat *, PetscReal **, PetscBool *a_isOK ); 257 #undef __FUNCT__ 258 #define __FUNCT__ "PCSetUp_GAMG" 259 PetscErrorCode PCSetUp_GAMG( PC pc ) 260 { 261 PetscErrorCode ierr; 262 PC_MG *mg = (PC_MG*)pc->data; 263 PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; 264 Mat Amat = pc->mat, Pmat = pc->pmat; 265 PetscBool isSeq, isMPI; 266 PetscInt fine_level, level, level1, M, N, bs, lidx; 267 MPI_Comm wcomm = ((PetscObject)pc)->comm; 268 PetscMPIInt mype,npe,nactivepe; 269 270 PetscFunctionBegin; 271 ierr = MPI_Comm_rank(wcomm,&mype);CHKERRQ(ierr); 272 ierr = MPI_Comm_size(wcomm,&npe);CHKERRQ(ierr); 273 if (pc->setupcalled){ 274 /* no state data in GAMG to destroy (now) */ 275 ierr = PCReset_MG(pc); CHKERRQ(ierr); 276 } 277 if (!pc_gamg->m_data) SETERRQ(wcomm,PETSC_ERR_SUP,"PCSetUp_GAMG called before PCSetCoordinates"); 278 /* setup special features of PCGAMG */ 279 ierr = PetscTypeCompare((PetscObject)Amat, MATSEQAIJ, &isSeq);CHKERRQ(ierr); 280 ierr = PetscTypeCompare((PetscObject)Amat, MATMPIAIJ, &isMPI);CHKERRQ(ierr); 281 if (isMPI) { 282 } else if (isSeq) { 283 } else SETERRQ1(wcomm,PETSC_ERR_ARG_WRONG, "Matrix type '%s' cannot be used with GAMG. GAMG can only handle AIJ matrices.",((PetscObject)Amat)->type_name); 284 285 /* GAMG requires input of fine-grid matrix. It determines nlevels. */ 286 ierr = MatGetBlockSize( Amat, &bs ); CHKERRQ(ierr); 287 if(bs!=1) SETERRQ1(wcomm,PETSC_ERR_ARG_WRONG, "GAMG only supports scalar prblems bs = '%d'.",bs); 288 289 /* Get A_i and R_i */ 290 ierr = MatGetSize( Amat, &M, &N );CHKERRQ(ierr); 291 PetscPrintf(PETSC_COMM_WORLD,"[%d]%s level %d N=%d\n",0,__FUNCT__,0,N); 292 #define GAMG_MAXLEVELS 20 293 Mat Aarr[GAMG_MAXLEVELS], Rarr[GAMG_MAXLEVELS]; PetscReal *coarse_crds = 0, *crds = pc_gamg->m_data; 294 PetscBool isOK; 295 for (level=0, Aarr[0] = Pmat, nactivepe = npe; level < GAMG_MAXLEVELS-1; level++ ){ 296 if( nactivepe == 1 ) { 297 break; 298 } 299 level1 = level + 1; 300 ierr = createProlongation( Aarr[level], crds, pc_gamg->m_dim, 301 &Rarr[level1], &coarse_crds, &isOK ); 302 CHKERRQ(ierr); 303 ierr = PetscFree( crds ); CHKERRQ( ierr ); 304 if(level==0) Aarr[0] = Amat; /* use Pmat for finest level setup, but use mat for solver */ 305 if( isOK ) { 306 ierr = partitionLevel( Aarr[level], pc_gamg->m_dim, &Rarr[level1], &coarse_crds, &nactivepe, &Aarr[level1] ); 307 CHKERRQ(ierr); 308 ierr = MatGetSize( Aarr[level1], &M, &N );CHKERRQ(ierr); 309 PetscPrintf(PETSC_COMM_WORLD,"[%d]%s done level %d N=%d, active pe = %d\n",0,__FUNCT__,level1,N,nactivepe); 310 } 311 else{ 312 break; 313 } 314 crds = coarse_crds; 315 } 316 ierr = PetscFree( coarse_crds ); CHKERRQ( ierr ); 317 PetscPrintf(PETSC_COMM_WORLD,"\t[%d]%s %d levels\n",0,__FUNCT__,level + 1); 318 pc_gamg->m_data = 0; /* destroyed coordinate data */ 319 pc_gamg->m_Nlevels = level + 1; 320 fine_level = level; 321 ierr = PCMGSetLevels(pc,pc_gamg->m_Nlevels,PETSC_NULL);CHKERRQ(ierr); 322 323 /* set default smoothers */ 324 PetscReal emax = 2.0, emin; 325 for (level=1; level<=fine_level; level++){ 326 KSP smoother; PC subpc; 327 ierr = PCMGGetSmoother(pc,level,&smoother);CHKERRQ(ierr); 328 ierr = KSPSetType(smoother,KSPCHEBYCHEV);CHKERRQ(ierr); 329 emin = emax/10.0; /* fix!!! */ 330 ierr = KSPGetPC(smoother,&subpc);CHKERRQ(ierr); 331 ierr = PCSetType(subpc,PCJACOBI);CHKERRQ(ierr); 332 ierr = KSPChebychevSetEigenvalues(smoother, emax, emin);CHKERRQ(ierr); /* need auto !!!!*/ 333 } 334 ierr = PCSetFromOptions_MG(pc); CHKERRQ(ierr); /* should be called in PCSetFromOptions_GAMG(), but cannot be called prior to PCMGSetLevels() */ 335 { 336 PetscBool galerkin; 337 ierr = PCMGGetGalerkin( pc, &galerkin); CHKERRQ(ierr); 338 if(galerkin){ 339 SETERRQ(wcomm,PETSC_ERR_ARG_WRONG, "GAMG does galerkin manually so it must not be used in PC_MG."); 340 } 341 } 342 /* create coarse level and the interpolation between the levels */ 343 for (level=0,lidx=pc_gamg->m_Nlevels-1; level<fine_level; level++,lidx--){ 344 level1 = level + 1; 345 ierr = PCMGSetInterpolation(pc,level1,Rarr[lidx]);CHKERRQ(ierr); 346 if(!PETSC_TRUE) { 347 PetscViewer viewer; char fname[32]; 348 sprintf(fname,"Amat_%d.m",lidx); 349 ierr = PetscViewerASCIIOpen( wcomm, fname, &viewer ); CHKERRQ(ierr); 350 ierr = PetscViewerSetFormat( viewer, PETSC_VIEWER_ASCII_MATLAB); CHKERRQ(ierr); 351 ierr = MatView( Aarr[lidx], viewer ); CHKERRQ(ierr); 352 ierr = PetscViewerDestroy( &viewer ); 353 } 354 ierr = MatDestroy( &Rarr[lidx] ); CHKERRQ(ierr); 355 { 356 KSP smoother; 357 ierr = PCMGGetSmoother(pc,level,&smoother); CHKERRQ(ierr); 358 ierr = KSPSetOperators( smoother, Aarr[lidx], Aarr[lidx], DIFFERENT_NONZERO_PATTERN ); 359 CHKERRQ(ierr); 360 ierr = MatDestroy( &Aarr[lidx] ); CHKERRQ(ierr); 361 } 362 } 363 { /* fine level (no P) */ 364 KSP smoother; 365 ierr = PCMGGetSmoother(pc,fine_level,&smoother); CHKERRQ(ierr); 366 ierr = KSPSetOperators( smoother, Aarr[0], Aarr[0], DIFFERENT_NONZERO_PATTERN ); 367 CHKERRQ(ierr); 368 } 369 370 /* setupcalled is set to 0 so that MG is setup from scratch */ 371 pc->setupcalled = 0; 372 ierr = PCSetUp_MG(pc);CHKERRQ(ierr); 373 PetscFunctionReturn(0); 374 } 375 376 /* -------------------------------------------------------------------------- */ 377 /* 378 PCDestroy_GAMG - Destroys the private context for the GAMG preconditioner 379 that was created with PCCreate_GAMG(). 380 381 Input Parameter: 382 . pc - the preconditioner context 383 384 Application Interface Routine: PCDestroy() 385 */ 386 #undef __FUNCT__ 387 #define __FUNCT__ "PCDestroy_GAMG" 388 PetscErrorCode PCDestroy_GAMG(PC pc) 389 { 390 PetscErrorCode ierr; 391 PC_MG *mg = (PC_MG*)pc->data; 392 PC_GAMG *pc_gamg= (PC_GAMG*)mg->innerctx; 393 394 PetscFunctionBegin; 395 ierr = PCReset_GAMG(pc);CHKERRQ(ierr); 396 ierr = PetscFree(pc_gamg);CHKERRQ(ierr); 397 ierr = PCDestroy_MG(pc);CHKERRQ(ierr); 398 PetscFunctionReturn(0); 399 } 400 401 #undef __FUNCT__ 402 #define __FUNCT__ "PCSetFromOptions_GAMG" 403 PetscErrorCode PCSetFromOptions_GAMG(PC pc) 404 { 405 /* PetscErrorCode ierr; */ 406 /* PC_MG *mg = (PC_MG*)pc->data; */ 407 /* PC_GAMG *pc_gamg = (PC_GAMG*)mg->innerctx; */ 408 /* MPI_Comm comm = ((PetscObject)pc)->comm; */ 409 410 PetscFunctionBegin; 411 PetscFunctionReturn(0); 412 } 413 414 /* -------------------------------------------------------------------------- */ 415 /* 416 PCCreate_GAMG - Creates a GAMG preconditioner context, PC_GAMG 417 418 Input Parameter: 419 . pc - the preconditioner context 420 421 Application Interface Routine: PCCreate() 422 423 */ 424 /* MC 425 PCGAMG - Use algebraic multigrid preconditioning. This preconditioner requires you provide 426 fine grid discretization matrix and coordinates on the fine grid. 427 428 Options Database Key: 429 Multigrid options(inherited) 430 + -pc_mg_cycles <1>: 1 for V cycle, 2 for W-cycle (MGSetCycles) 431 . -pc_mg_smoothup <1>: Number of post-smoothing steps (MGSetNumberSmoothUp) 432 . -pc_mg_smoothdown <1>: Number of pre-smoothing steps (MGSetNumberSmoothDown) 433 -pc_mg_type <multiplicative>: (one of) additive multiplicative full cascade kascade 434 GAMG options: 435 436 Level: intermediate 437 Concepts: multigrid 438 439 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, PCMGType, 440 PCMGSetLevels(), PCMGGetLevels(), PCMGSetType(), MPSetCycles(), PCMGSetNumberSmoothDown(), 441 PCMGSetNumberSmoothUp(), PCMGGetCoarseSolve(), PCMGSetResidual(), PCMGSetInterpolation(), 442 PCMGSetRestriction(), PCMGGetSmoother(), PCMGGetSmootherUp(), PCMGGetSmootherDown(), 443 PCMGSetCyclesOnLevel(), PCMGSetRhs(), PCMGSetX(), PCMGSetR() 444 M */ 445 446 EXTERN_C_BEGIN 447 #undef __FUNCT__ 448 #define __FUNCT__ "PCCreate_GAMG" 449 PetscErrorCode PCCreate_GAMG(PC pc) 450 { 451 PetscErrorCode ierr; 452 PC_GAMG *pc_gamg; 453 PC_MG *mg; 454 455 PetscFunctionBegin; 456 /* PCGAMG is an inherited class of PCMG. Initialize pc as PCMG */ 457 ierr = PCSetType(pc,PCMG);CHKERRQ(ierr); /* calls PCCreate_MG() and MGCreate_Private() */ 458 ierr = PetscObjectChangeTypeName((PetscObject)pc,PCGAMG);CHKERRQ(ierr); 459 460 /* create a supporting struct and attach it to pc */ 461 ierr = PetscNewLog(pc,PC_GAMG,&pc_gamg);CHKERRQ(ierr); 462 mg = (PC_MG*)pc->data; 463 mg->innerctx = pc_gamg; 464 465 pc_gamg->m_Nlevels = -1; 466 467 /* overwrite the pointers of PCMG by the functions of PCGAMG */ 468 pc->ops->setfromoptions = PCSetFromOptions_GAMG; 469 pc->ops->setup = PCSetUp_GAMG; 470 pc->ops->reset = PCReset_GAMG; 471 pc->ops->destroy = PCDestroy_GAMG; 472 473 ierr = PetscObjectComposeFunctionDynamic( (PetscObject)pc, 474 "PCSetCoordinates_C", 475 "PCSetCoordinates_GAMG", 476 PCSetCoordinates_GAMG);CHKERRQ(ierr); 477 PetscFunctionReturn(0); 478 } 479 EXTERN_C_END 480