Actual source code: damg.c

      1: #define PETSCSNES_DLL
      2:
     3:  #include petscda.h
     4:  #include petscksp.h
     5:  #include petscmg.h
     6:  #include petscdmmg.h
     7:  #include src/ksp/pc/pcimpl.h

      9: /*
     10:    Code for almost fully managing multigrid/multi-level linear solvers for DA grids
     11: */

     15: /*@C
     16:     DMMGCreate - Creates a DA based multigrid solver object. This allows one to
     17:       easily implement MG methods on regular grids.

     19:     Collective on MPI_Comm

     21:     Input Parameter:
     22: +   comm - the processors that will share the grids and solution process
     23: .   nlevels - number of multigrid levels
     24: -   user - an optional user context

     26:     Output Parameters:
     27: .    - the context

     29:     Notes:
     30:       To provide a different user context for each level call DMMGSetUser() after calling
     31:       this routine

     33:     Level: advanced

     35: .seealso DMMGDestroy(), DMMGSetUser(), DMMGGetUser()

     37: @*/
     38: PetscErrorCode PETSCSNES_DLLEXPORT DMMGCreate(MPI_Comm comm,PetscInt nlevels,void *user,DMMG **dmmg)
     39: {
     41:   PetscInt       i;
     42:   DMMG           *p;
     43:   PetscTruth     galerkin;

     46:   PetscOptionsGetInt(0,"-dmmg_nlevels",&nlevels,PETSC_IGNORE);
     47:   PetscOptionsHasName(0,"-dmmg_galerkin",&galerkin);

     49:   PetscMalloc(nlevels*sizeof(DMMG),&p);
     50:   for (i=0; i<nlevels; i++) {
     51:     PetscNew(struct _p_DMMG,&p[i]);
     52:     p[i]->nlevels  = nlevels - i;
     53:     p[i]->comm     = comm;
     54:     p[i]->user     = user;
     55:     p[i]->galerkin = galerkin;
     56:   }
     57:   *dmmg = p;
     58:   return(0);
     59: }

     63: /*@C
     64:     DMMGSetUseGalerkinCoarse - Courses the DMMG to use R*A_f*R^T to form
     65:        the coarser matrices from finest

     67:     Collective on DMMG

     69:     Input Parameter:
     70: .    - the context

     72:     Options Database Keys:
     73: .    -dmmg_galerkin

     75:     Level: advanced

     77: .seealso DMMGCreate()

     79: @*/
     80: PetscErrorCode PETSCSNES_DLLEXPORT DMMGSetUseGalerkinCoarse(DMMG* dmmg)
     81: {
     82:   PetscInt  i,nlevels = dmmg[0]->nlevels;

     85:   if (!dmmg) SETERRQ(PETSC_ERR_ARG_NULL,"Passing null as DMMG");

     87:   for (i=0; i<nlevels; i++) {
     88:     dmmg[i]->galerkin = PETSC_TRUE;
     89:   }
     90:   return(0);
     91: }

     95: /*@C
     96:     DMMGDestroy - Destroys a DA based multigrid solver object.

     98:     Collective on DMMG

    100:     Input Parameter:
    101: .    - the context

    103:     Level: advanced

    105: .seealso DMMGCreate()

    107: @*/
    108: PetscErrorCode PETSCSNES_DLLEXPORT DMMGDestroy(DMMG *dmmg)
    109: {
    111:   PetscInt       i,nlevels = dmmg[0]->nlevels;

    114:   if (!dmmg) SETERRQ(PETSC_ERR_ARG_NULL,"Passing null as DMMG");

    116:   for (i=1; i<nlevels; i++) {
    117:     if (dmmg[i]->R) {MatDestroy(dmmg[i]->R);}
    118:   }
    119:   for (i=0; i<nlevels; i++) {
    120:     if (dmmg[i]->dm)      {DMDestroy(dmmg[i]->dm);}
    121:     if (dmmg[i]->x)       {VecDestroy(dmmg[i]->x);}
    122:     if (dmmg[i]->b)       {VecDestroy(dmmg[i]->b);}
    123:     if (dmmg[i]->r)       {VecDestroy(dmmg[i]->r);}
    124:     if (dmmg[i]->work1)   {VecDestroy(dmmg[i]->work1);}
    125:     if (dmmg[i]->w)       {VecDestroy(dmmg[i]->w);}
    126:     if (dmmg[i]->work2)   {VecDestroy(dmmg[i]->work2);}
    127:     if (dmmg[i]->lwork1)  {VecDestroy(dmmg[i]->lwork1);}
    128:     if (dmmg[i]->B && dmmg[i]->B != dmmg[i]->J) {MatDestroy(dmmg[i]->B);}
    129:     if (dmmg[i]->J)         {MatDestroy(dmmg[i]->J);}
    130:     if (dmmg[i]->Rscale)    {VecDestroy(dmmg[i]->Rscale);}
    131:     if (dmmg[i]->fdcoloring){MatFDColoringDestroy(dmmg[i]->fdcoloring);}
    132:     if (dmmg[i]->ksp && !dmmg[i]->snes) {KSPDestroy(dmmg[i]->ksp);}
    133:     if (dmmg[i]->snes)      {PetscObjectDestroy((PetscObject)dmmg[i]->snes);}
    134:     if (dmmg[i]->inject)    {VecScatterDestroy(dmmg[i]->inject);}
    135:     PetscFree(dmmg[i]);
    136:   }
    137:   PetscFree(dmmg);
    138:   return(0);
    139: }

    143: /*@C
    144:     DMMGSetDM - Sets the coarse grid information for the grids

    146:     Collective on DMMG

    148:     Input Parameter:
    149: +   dmmg - the context
    150: -   dm - the DA or VecPack object

    152:     Level: advanced

    154: .seealso DMMGCreate(), DMMGDestroy()

    156: @*/
    157: PetscErrorCode PETSCSNES_DLLEXPORT DMMGSetDM(DMMG *dmmg,DM dm)
    158: {
    160:   PetscInt       i,nlevels = dmmg[0]->nlevels;

    163:   if (!dmmg) SETERRQ(PETSC_ERR_ARG_NULL,"Passing null as DMMG");

    165:   /* Create DA data structure for all the levels */
    166:   dmmg[0]->dm = dm;
    167:   PetscObjectReference((PetscObject)dm);
    168:   for (i=1; i<nlevels; i++) {
    169:     DMRefine(dmmg[i-1]->dm,dmmg[i]->comm,&dmmg[i]->dm);
    170:   }
    171:   DMMGSetUp(dmmg);
    172:   return(0);
    173: }

    177: /*@C
    178:     DMMGSetUp - Prepares the DMMG to solve a system

    180:     Collective on DMMG

    182:     Input Parameter:
    183: .   dmmg - the context

    185:     Level: advanced

    187: .seealso DMMGCreate(), DMMGDestroy(), DMMG, DMMGSetSNES(), DMMGSetKSP(), DMMGSolve()

    189: @*/
    190: PetscErrorCode PETSCSNES_DLLEXPORT DMMGSetUp(DMMG *dmmg)
    191: {
    193:   PetscInt       i,nlevels = dmmg[0]->nlevels;

    197:   /* Create work vectors and matrix for each level */
    198:   for (i=0; i<nlevels; i++) {
    199:     DMCreateGlobalVector(dmmg[i]->dm,&dmmg[i]->x);
    200:     VecDuplicate(dmmg[i]->x,&dmmg[i]->b);
    201:     VecDuplicate(dmmg[i]->x,&dmmg[i]->r);
    202:   }

    204:   /* Create interpolation/restriction between levels */
    205:   for (i=1; i<nlevels; i++) {
    206:     DMGetInterpolation(dmmg[i-1]->dm,dmmg[i]->dm,&dmmg[i]->R,PETSC_NULL);
    207:   }

    209:   return(0);
    210: }

    214: /*@C
    215:     DMMGSolve - Actually solves the (non)linear system defined with the DMMG

    217:     Collective on DMMG

    219:     Input Parameter:
    220: .   dmmg - the context

    222:     Level: advanced

    224:     Options Database:
    225: +   -dmmg_grid_sequence - use grid sequencing to get the initial solution for each level from the previous
    226: -   -dmmg_vecmonitor - display the solution at each iteration

    228:      Notes: For linear (KSP) problems may be called more than once, uses the same
    229:     matrices but recomputes the right hand side for each new solve. Call DMMGSetKSP()
    230:     to generate new matrices.
    231:
    232: .seealso DMMGCreate(), DMMGDestroy(), DMMG, DMMGSetSNES(), DMMGSetKSP(), DMMGSetUp()

    234: @*/
    235: PetscErrorCode PETSCSNES_DLLEXPORT DMMGSolve(DMMG *dmmg)
    236: {
    238:   PetscInt       i,nlevels = dmmg[0]->nlevels;
    239:   PetscTruth     gridseq,vecmonitor,flg;

    242:   PetscOptionsHasName(0,"-dmmg_grid_sequence",&gridseq);
    243:   PetscOptionsHasName(0,"-dmmg_vecmonitor",&vecmonitor);
    244:   if (gridseq) {
    245:     if (dmmg[0]->initialguess) {
    246:       (*dmmg[0]->initialguess)(dmmg[0],dmmg[0]->x);
    247:       if (dmmg[0]->ksp && !dmmg[0]->snes) {
    248:         KSPSetInitialGuessNonzero(dmmg[0]->ksp,PETSC_TRUE);
    249:       }
    250:     }
    251:     for (i=0; i<nlevels-1; i++) {
    252:       (*dmmg[i]->solve)(dmmg,i);
    253:       if (vecmonitor) {
    254:         VecView(dmmg[i]->x,PETSC_VIEWER_DRAW_(dmmg[i]->comm));
    255:       }
    256:       MatInterpolate(dmmg[i+1]->R,dmmg[i]->x,dmmg[i+1]->x);
    257:       if (dmmg[i+1]->ksp && !dmmg[i+1]->ksp) {
    258:         KSPSetInitialGuessNonzero(dmmg[i+1]->ksp,PETSC_TRUE);
    259:       }
    260:     }
    261:   } else {
    262:     if (dmmg[nlevels-1]->initialguess) {
    263:       (*dmmg[nlevels-1]->initialguess)(dmmg[nlevels-1],dmmg[nlevels-1]->x);
    264:     }
    265:   }
    266:   (*DMMGGetFine(dmmg)->solve)(dmmg,nlevels-1);
    267:   if (vecmonitor) {
    268:      VecView(dmmg[nlevels-1]->x,PETSC_VIEWER_DRAW_(dmmg[nlevels-1]->comm));
    269:   }

    271:   PetscOptionsHasName(PETSC_NULL,"-dmmg_view",&flg);
    272:   if (flg && !PetscPreLoadingOn) {
    273:     DMMGView(dmmg,PETSC_VIEWER_STDOUT_(dmmg[0]->comm));
    274:   }
    275:   PetscOptionsHasName(PETSC_NULL,"-dmmg_view_binary",&flg);
    276:   if (flg && !PetscPreLoadingOn) {
    277:     DMMGView(dmmg,PETSC_VIEWER_BINARY_(dmmg[0]->comm));
    278:   }
    279:   return(0);
    280: }

    284: PetscErrorCode PETSCSNES_DLLEXPORT DMMGSolveKSP(DMMG *dmmg,PetscInt level)
    285: {

    289:   if (dmmg[level]->rhs) {
    290:     (*dmmg[level]->rhs)(dmmg[level],dmmg[level]->b);
    291:   }
    292:   if (dmmg[level]->matricesset) {
    293:     KSPSetOperators(dmmg[level]->ksp,dmmg[level]->J,dmmg[level]->J,SAME_NONZERO_PATTERN);
    294:     dmmg[level]->matricesset = PETSC_FALSE;
    295:   }
    296:   KSPSolve(dmmg[level]->ksp,dmmg[level]->b,dmmg[level]->x);
    297:   return(0);
    298: }

    300: /*
    301:     Sets each of the linear solvers to use multigrid
    302: */
    305: PetscErrorCode PETSCSNES_DLLEXPORT DMMGSetUpLevel(DMMG *dmmg,KSP ksp,PetscInt nlevels)
    306: {
    308:   PetscInt       i;
    309:   PC             pc;
    310:   PetscTruth     ismg,monitor,ismf,isshell,ismffd;
    311:   KSP            lksp; /* solver internal to the multigrid preconditioner */
    312:   MPI_Comm       *comms,comm;
    313:   PetscViewer    ascii;

    316:   if (!dmmg) SETERRQ(PETSC_ERR_ARG_NULL,"Passing null as DMMG");

    318:   PetscOptionsHasName(PETSC_NULL,"-dmmg_ksp_monitor",&monitor);
    319:   if (monitor) {
    320:     PetscObjectGetComm((PetscObject)ksp,&comm);
    321:     PetscViewerASCIIOpen(comm,"stdout",&ascii);
    322:     PetscViewerASCIISetTab(ascii,1+dmmg[0]->nlevels-nlevels);
    323:     KSPSetMonitor(ksp,KSPDefaultMonitor,ascii,(PetscErrorCode(*)(void*))PetscViewerDestroy);
    324:   }

    326:   /* use fgmres on outer iteration by default */
    327:   KSPSetType(ksp,KSPFGMRES);
    328:   KSPGetPC(ksp,&pc);
    329:   PCSetType(pc,PCMG);
    330:   PetscMalloc(nlevels*sizeof(MPI_Comm),&comms);
    331:   for (i=0; i<nlevels; i++) {
    332:     comms[i] = dmmg[i]->comm;
    333:   }
    334:   PCMGSetLevels(pc,nlevels,comms);
    335:   PetscFree(comms);
    336:    PCMGSetType(pc,PC_MG_FULL);

    338:   PetscTypeCompare((PetscObject)pc,PCMG,&ismg);
    339:   if (ismg) {

    341:     /* set solvers for each level */
    342:     for (i=0; i<nlevels; i++) {
    343:       PCMGGetSmoother(pc,i,&lksp);
    344:       KSPSetOperators(lksp,dmmg[i]->J,dmmg[i]->B,DIFFERENT_NONZERO_PATTERN);
    345:       if (i < nlevels-1) { /* don't set for finest level, they are set in PCApply_MG()*/
    346:         PCMGSetX(pc,i,dmmg[i]->x);
    347:         PCMGSetRhs(pc,i,dmmg[i]->b);
    348:       }
    349:       if (i > 0) {
    350:         PCMGSetR(pc,i,dmmg[i]->r);
    351:         PCMGSetResidual(pc,i,PCMGDefaultResidual,dmmg[i]->J);
    352:       }
    353:       if (monitor) {
    354:         PetscObjectGetComm((PetscObject)lksp,&comm);
    355:         PetscViewerASCIIOpen(comm,"stdout",&ascii);
    356:         PetscViewerASCIISetTab(ascii,1+dmmg[0]->nlevels-i);
    357:         KSPSetMonitor(lksp,KSPDefaultMonitor,ascii,(PetscErrorCode(*)(void*))PetscViewerDestroy);
    358:       }
    359:       /* If using a matrix free multiply and did not provide an explicit matrix to build
    360:          the preconditioner then must use no preconditioner
    361:       */
    362:       PetscTypeCompare((PetscObject)dmmg[i]->B,MATSHELL,&isshell);
    363:       PetscTypeCompare((PetscObject)dmmg[i]->B,MATDAAD,&ismf);
    364:       PetscTypeCompare((PetscObject)dmmg[i]->B,MATMFFD,&ismffd);
    365:       if (isshell || ismf || ismffd) {
    366:         PC  lpc;
    367:         KSPGetPC(lksp,&lpc);
    368:         PCSetType(lpc,PCNONE);
    369:       }
    370:     }

    372:     /* Set interpolation/restriction between levels */
    373:     for (i=1; i<nlevels; i++) {
    374:       PCMGSetInterpolate(pc,i,dmmg[i]->R);
    375:       PCMGSetRestriction(pc,i,dmmg[i]->R);
    376:     }
    377:   }
    378:   return(0);
    379: }

    383: /*@C
    384:     DMMGSetKSP - Sets the linear solver object that will use the grid hierarchy

    386:     Collective on DMMG

    388:     Input Parameter:
    389: +   dmmg - the context
    390: .   func - function to compute linear system matrix on each grid level
    391: -   rhs - function to compute right hand side on each level (need only work on the finest grid
    392:           if you do not use grid sequencing)

    394:     Level: advanced

    396:     Notes: For linear problems my be called more than once, reevaluates the matrices if it is called more
    397:        than once. Call DMMGSolve() directly several times to solve with the same matrix but different
    398:        right hand sides.
    399:
    400: .seealso DMMGCreate(), DMMGDestroy, DMMGSetDM(), DMMGSolve()

    402: @*/
    403: PetscErrorCode PETSCSNES_DLLEXPORT DMMGSetKSP(DMMG *dmmg,PetscErrorCode (*rhs)(DMMG,Vec),PetscErrorCode (*func)(DMMG,Mat))
    404: {
    406:   PetscInt       i,nlevels = dmmg[0]->nlevels;
    407:   PetscTruth     galerkin;

    410:   if (!dmmg) SETERRQ(PETSC_ERR_ARG_NULL,"Passing null as DMMG");
    411:   galerkin = dmmg[0]->galerkin;

    413:   if (galerkin) {
    414:     DMGetMatrix(dmmg[nlevels-1]->dm,MATAIJ,&dmmg[nlevels-1]->B);
    415:     (*func)(dmmg[nlevels-1],dmmg[nlevels-1]->B);
    416:     for (i=nlevels-2; i>-1; i--) {
    417:       MatPtAP(dmmg[i+1]->B,dmmg[i+1]->R,MAT_INITIAL_MATRIX,1.0,&dmmg[i]->B);
    418:     }
    419:   }

    421:   if (!dmmg[0]->ksp) {
    422:     /* create solvers for each level if they don't already exist*/
    423:     for (i=0; i<nlevels; i++) {

    425:       if (!dmmg[i]->B && !galerkin) {
    426:         DMGetMatrix(dmmg[i]->dm,MATAIJ,&dmmg[i]->B);
    427:       }
    428:       if (!dmmg[i]->J) {
    429:         dmmg[i]->J = dmmg[i]->B;
    430:       }

    432:       KSPCreate(dmmg[i]->comm,&dmmg[i]->ksp);
    433:       DMMGSetUpLevel(dmmg,dmmg[i]->ksp,i+1);
    434:       KSPSetFromOptions(dmmg[i]->ksp);
    435:       dmmg[i]->solve = DMMGSolveKSP;
    436:       dmmg[i]->rhs   = rhs;
    437:     }
    438:   }

    440:   /* evalute matrix on each level */
    441:   for (i=0; i<nlevels; i++) {
    442:     if (!galerkin) {
    443:       (*func)(dmmg[i],dmmg[i]->J);
    444:     }
    445:     dmmg[i]->matricesset = PETSC_TRUE;
    446:   }

    448:   for (i=0; i<nlevels-1; i++) {
    449:     KSPSetOptionsPrefix(dmmg[i]->ksp,"dmmg_");
    450:   }

    452:   return(0);
    453: }

    457: /*@C
    458:     DMMGView - prints information on a DA based multi-level preconditioner

    460:     Collective on DMMG and PetscViewer

    462:     Input Parameter:
    463: +   dmmg - the context
    464: -   viewer - the viewer

    466:     Level: advanced

    468: .seealso DMMGCreate(), DMMGDestroy

    470: @*/
    471: PetscErrorCode PETSCSNES_DLLEXPORT DMMGView(DMMG *dmmg,PetscViewer viewer)
    472: {
    474:   PetscInt       i,nlevels = dmmg[0]->nlevels;
    475:   PetscMPIInt    flag;
    476:   MPI_Comm       comm;
    477:   PetscTruth     iascii,isbinary;

    482:   PetscObjectGetComm((PetscObject)viewer,&comm);
    483:   MPI_Comm_compare(comm,dmmg[0]->comm,&flag);
    484:   if (flag != MPI_CONGRUENT && flag != MPI_IDENT) {
    485:     SETERRQ(PETSC_ERR_ARG_NOTSAMECOMM,"Different communicators in the DMMG and the PetscViewer");
    486:   }

    488:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
    489:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_BINARY,&isbinary);
    490:   if (isbinary) {
    491:     for (i=0; i<nlevels; i++) {
    492:       MatView(dmmg[i]->J,viewer);
    493:     }
    494:     for (i=1; i<nlevels; i++) {
    495:       MatView(dmmg[i]->R,viewer);
    496:     }
    497:   } else {
    498:     if (iascii) {
    499:       PetscViewerASCIIPrintf(viewer,"DMMG Object with %D levels\n",nlevels);
    500:     }
    501:     for (i=0; i<nlevels; i++) {
    502:       PetscViewerASCIIPushTab(viewer);
    503:       DMView(dmmg[i]->dm,viewer);
    504:       PetscViewerASCIIPopTab(viewer);
    505:     }
    506:     if (iascii) {
    507:       PetscViewerASCIIPrintf(viewer,"%s Object on finest level\n",dmmg[nlevels-1]->ksp ? "KSP" : "SNES");
    508:       if (dmmg[nlevels-1]->galerkin) {
    509:         PetscViewerASCIIPrintf(viewer,"Using Galerkin R^T*A*R process to compute coarser matrices");
    510:       }
    511:     }
    512:     if (dmmg[nlevels-1]->ksp) {
    513:       KSPView(dmmg[nlevels-1]->ksp,viewer);
    514:     } else {
    515:       /* use of PetscObjectView() means we do not have to link with libpetscsnes if SNES is not being used */
    516:       PetscObjectView((PetscObject)dmmg[nlevels-1]->snes,viewer);
    517:     }
    518:   }
    519:   return(0);
    520: }

    524: /*@C
    525:     DMMGSetNullSpace - Indicates the null space in the linear operator (this is needed by the linear solver)

    527:     Collective on DMMG

    529:     Input Parameter:
    530: +   dmmg - the context
    531: .   has_cnst - is the constant vector in the null space
    532: .   n - number of null vectors (excluding the possible constant vector)
    533: -   func - a function that fills an array of vectors with the null vectors (must be orthonormal), may be PETSC_NULL

    535:     Level: advanced

    537: .seealso DMMGCreate(), DMMGDestroy, DMMGSetDM(), DMMGSolve(), MatNullSpaceCreate(), KSPSetNullSpace()

    539: @*/
    540: PetscErrorCode PETSCSNES_DLLEXPORT DMMGSetNullSpace(DMMG *dmmg,PetscTruth has_cnst,PetscInt n,PetscErrorCode (*func)(DMMG,Vec[]))
    541: {
    543:   PetscInt       i,j,nlevels = dmmg[0]->nlevels;
    544:   Vec            *nulls = 0;
    545:   MatNullSpace   nullsp;
    546:   KSP            iksp;
    547:   PC             pc,ipc;
    548:   PetscTruth     ismg,isred;

    551:   if (!dmmg) SETERRQ(PETSC_ERR_ARG_NULL,"Passing null as DMMG");
    552:   if (!dmmg[0]->ksp) SETERRQ(PETSC_ERR_ORDER,"Must call AFTER DMMGSetKSP() or DMMGSetSNES()");
    553:   if ((n && !func) || (!n && func)) SETERRQ(PETSC_ERR_ARG_INCOMP,"Both n and func() must be set together");
    554:   if (n < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Cannot have negative number of vectors in null space n = %D",n)

    556:   for (i=0; i<nlevels; i++) {
    557:     if (n) {
    558:       VecDuplicateVecs(dmmg[i]->b,n,&nulls);
    559:       (*func)(dmmg[i],nulls);
    560:     }
    561:     MatNullSpaceCreate(dmmg[i]->comm,has_cnst,n,nulls,&nullsp);
    562:     KSPSetNullSpace(dmmg[i]->ksp,nullsp);
    563:     for (j=i; j<nlevels; j++) {
    564:       KSPGetPC(dmmg[j]->ksp,&pc);
    565:       PetscTypeCompare((PetscObject)pc,PCMG,&ismg);
    566:       if (ismg) {
    567:         PCMGGetSmoother(pc,i,&iksp);
    568:         KSPSetNullSpace(iksp, nullsp);
    569:       }
    570:     }
    571:     MatNullSpaceDestroy(nullsp);
    572:     if (n) {
    573:       PetscFree(nulls);
    574:     }
    575:   }
    576:   /* make all the coarse grid solvers have LU shift since they are singular */
    577:   for (i=0; i<nlevels; i++) {
    578:     KSPGetPC(dmmg[i]->ksp,&pc);
    579:     PetscTypeCompare((PetscObject)pc,PCMG,&ismg);
    580:     if (ismg) {
    581:       PCMGGetSmoother(pc,0,&iksp);
    582:       KSPGetPC(iksp,&ipc);
    583:       PetscTypeCompare((PetscObject)ipc,PCREDUNDANT,&isred);
    584:       if (isred) {
    585:         PCRedundantGetPC(ipc,&ipc);
    586:       }
    587:       PCFactorSetShiftPd(ipc,PETSC_TRUE);
    588:     }
    589:   }
    590:   return(0);
    591: }

    595: /*@C
    596:     DMMGInitialGuessCurrent - Use with DMMGSetInitialGuess() to use the current value in the
    597:        solution vector (obtainable with DMMGGetx() as the initial guess)

    599:     Collective on DMMG

    601:     Input Parameter:
    602: +   dmmg - the context
    603: -   vec - dummy argument

    605:     Level: intermediate

    607: .seealso DMMGCreate(), DMMGDestroy, DMMGSetKSP(), DMMGSetSNES(), DMMGSetInitialGuess()

    609: @*/
    610: PetscErrorCode PETSCSNES_DLLEXPORT DMMGInitialGuessCurrent(DMMG dmmg,Vec vec)
    611: {
    613:   return(0);
    614: }

    618: /*@C
    619:     DMMGSetInitialGuess - Sets the function that computes an initial guess.

    621:     Collective on DMMG

    623:     Input Parameter:
    624: +   dmmg - the context
    625: -   guess - the function

    627:     Notes: For nonlinear problems, if this is not set, then the current value in the
    628:              solution vector (obtained with DMMGGetX()) is used. Thus is if you doing 'time
    629:              stepping' it will use your current solution as the guess for the next timestep.
    630:            If grid sequencing is used (via -dmmg_grid_sequence) then the "guess" function
    631:              is used only on the coarsest grid.
    632:            For linear problems, if this is not set, then 0 is used as an initial guess.
    633:              If you would like the linear solver to also (like the nonlinear solver) use
    634:              the current solution vector as the initial guess then use DMMGInitialGuessCurrent()
    635:              as the function you pass in

    637:     Level: intermediate

    640: .seealso DMMGCreate(), DMMGDestroy, DMMGSetKSP(), DMMGSetSNES(), DMMGInitialGuessCurrent()

    642: @*/
    643: PetscErrorCode PETSCSNES_DLLEXPORT DMMGSetInitialGuess(DMMG *dmmg,PetscErrorCode (*guess)(DMMG,Vec))
    644: {
    645:   PetscInt i,nlevels = dmmg[0]->nlevels;

    648:   for (i=0; i<nlevels; i++) {
    649:     dmmg[i]->initialguess = guess;
    650:   }
    651:   return(0);
    652: }