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