Actual source code: damg.c

  1: 
 2:  #include petscda.h
 3:  #include petscksp.h
 4:  #include petscmg.h

  6: /*
  7:    Code for almost fully managing multigrid/multi-level linear solvers for DA grids
  8: */

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

 16:     Collective on MPI_Comm

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

 23:     Output Parameters:
 24: .    - the context

 26:     Notes:
 27:       To provide a different user context for each level call DMMGSetUser() after calling
 28:       this routine

 30:     Level: advanced

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

 34: @*/
 35: PetscErrorCode DMMGCreate(MPI_Comm comm,PetscInt nlevels,void *user,DMMG **dmmg)
 36: {
 38:   PetscInt       i;
 39:   DMMG           *p;
 40:   PetscTruth     galerkin;

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

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

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

 65:     Collective on DMMG

 67:     Input Parameter:
 68: .    - the context

 70:     Options Database Keys:
 71: .    -dmmg_galerkin

 73:     Level: advanced

 75: .seealso DMMGCreate()

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

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

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

 93: /*@C
 94:     DMMGDestroy - Destroys a DA based multigrid solver object. 

 96:     Collective on DMMG

 98:     Input Parameter:
 99: .    - the context

101:     Level: advanced

103: .seealso DMMGCreate()

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

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

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

141: /*@C
142:     DMMGSetDM - Sets the coarse grid information for the grids

144:     Collective on DMMG

146:     Input Parameter:
147: +   dmmg - the context
148: -   dm - the DA or VecPack object

150:     Level: advanced

152: .seealso DMMGCreate(), DMMGDestroy()

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

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

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

175: /*@C
176:     DMMGSetUp - Prepares the DMMG to solve a system

178:     Collective on DMMG

180:     Input Parameter:
181: .   dmmg - the context

183:     Level: advanced

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

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


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

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

207:   return(0);
208: }

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

215:     Collective on DMMG

217:     Input Parameter:
218: .   dmmg - the context

220:     Level: advanced

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

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

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

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

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

282: PetscErrorCode DMMGSolveKSP(DMMG *dmmg,PetscInt level)
283: {

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

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

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

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

322:   /* use fgmres on outer iteration by default */
323:   KSPSetType(ksp,KSPFGMRES);
324:   KSPGetPC(ksp,&pc);
325:   PCSetType(pc,PCMG);
326:   PetscMalloc(nlevels*sizeof(MPI_Comm),&comms);
327:   for (i=0; i<nlevels; i++) {
328:     comms[i] = dmmg[i]->comm;
329:   }
330:   MGSetLevels(pc,nlevels,comms);
331:   PetscFree(comms);
332:    MGSetType(pc,MGFULL);

334:   PetscTypeCompare((PetscObject)pc,PCMG,&ismg);
335:   if (ismg) {

337:     /* set solvers for each level */
338:     for (i=0; i<nlevels; i++) {
339:       MGGetSmoother(pc,i,&lksp);
340:       KSPSetOperators(lksp,dmmg[i]->J,dmmg[i]->B,DIFFERENT_NONZERO_PATTERN);
341:       MGSetX(pc,i,dmmg[i]->x);
342:       MGSetRhs(pc,i,dmmg[i]->b);
343:       MGSetR(pc,i,dmmg[i]->r);
344:       MGSetResidual(pc,i,MGDefaultResidual,dmmg[i]->J);
345:       if (monitor) {
346:         PetscObjectGetComm((PetscObject)lksp,&comm);
347:         PetscViewerASCIIOpen(comm,"stdout",&ascii);
348:         PetscViewerASCIISetTab(ascii,1+dmmg[0]->nlevels-i);
349:         KSPSetMonitor(lksp,KSPDefaultMonitor,ascii,(PetscErrorCode(*)(void*))PetscViewerDestroy);
350:       }
351:       /* If using a matrix free multiply and did not provide an explicit matrix to build
352:          the preconditioner then must use no preconditioner 
353:       */
354:       PetscTypeCompare((PetscObject)dmmg[i]->B,MATSHELL,&isshell);
355:       PetscTypeCompare((PetscObject)dmmg[i]->B,MATDAAD,&ismf);
356:       PetscTypeCompare((PetscObject)dmmg[i]->B,MATMFFD,&ismffd);
357:       if (isshell || ismf || ismffd) {
358:         PC  lpc;
359:         KSPGetPC(lksp,&lpc);
360:         PCSetType(lpc,PCNONE);
361:       }
362:     }

364:     /* Set interpolation/restriction between levels */
365:     for (i=1; i<nlevels; i++) {
366:       MGSetInterpolate(pc,i,dmmg[i]->R);
367:       MGSetRestriction(pc,i,dmmg[i]->R);
368:     }
369:   }
370:   return(0);
371: }

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

378:     Collective on DMMG

380:     Input Parameter:
381: +   dmmg - the context
382: .   func - function to compute linear system matrix on each grid level
383: -   rhs - function to compute right hand side on each level (need only work on the finest grid
384:           if you do not use grid sequencing)

386:     Level: advanced

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

394: @*/
395: PetscErrorCode DMMGSetKSP(DMMG *dmmg,PetscErrorCode (*rhs)(DMMG,Vec),PetscErrorCode (*func)(DMMG,Mat))
396: {
398:   PetscInt       i,nlevels = dmmg[0]->nlevels;
399:   PetscTruth     galerkin;

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

405:   if (galerkin) {
406:     DMGetMatrix(dmmg[nlevels-1]->dm,MATAIJ,&dmmg[nlevels-1]->B);
407:     (*func)(dmmg[nlevels-1],dmmg[nlevels-1]->B);
408:     for (i=nlevels-2; i>-1; i--) {
409:       MatPtAP(dmmg[i+1]->B,dmmg[i+1]->R,MAT_INITIAL_MATRIX,1.0,&dmmg[i]->B);
410:     }
411:   }

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

417:       if (!dmmg[i]->B && !galerkin) {
418:         DMGetMatrix(dmmg[i]->dm,MATAIJ,&dmmg[i]->B);
419:       }
420:       if (!dmmg[i]->J) {
421:         dmmg[i]->J = dmmg[i]->B;
422:       }

424:       KSPCreate(dmmg[i]->comm,&dmmg[i]->ksp);
425:       DMMGSetUpLevel(dmmg,dmmg[i]->ksp,i+1);
426:       KSPSetFromOptions(dmmg[i]->ksp);
427:       dmmg[i]->solve = DMMGSolveKSP;
428:       dmmg[i]->rhs   = rhs;
429:     }
430:   }

432:   /* evalute matrix on each level */
433:   for (i=0; i<nlevels; i++) {
434:     if (!galerkin) {
435:       (*func)(dmmg[i],dmmg[i]->J);
436:     }
437:     dmmg[i]->matricesset = PETSC_TRUE;
438:   }

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

444:   return(0);
445: }

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

452:     Collective on DMMG and PetscViewer

454:     Input Parameter:
455: +   dmmg - the context
456: -   viewer - the viewer

458:     Level: advanced

460: .seealso DMMGCreate(), DMMGDestroy

462: @*/
463: PetscErrorCode DMMGView(DMMG *dmmg,PetscViewer viewer)
464: {
466:   PetscInt       i,nlevels = dmmg[0]->nlevels;
467:   PetscMPIInt    flag;
468:   MPI_Comm       comm;
469:   PetscTruth     iascii,isbinary;

474:   PetscObjectGetComm((PetscObject)viewer,&comm);
475:   MPI_Comm_compare(comm,dmmg[0]->comm,&flag);
476:   if (flag != MPI_CONGRUENT && flag != MPI_IDENT) {
477:     SETERRQ(PETSC_ERR_ARG_NOTSAMECOMM,"Different communicators in the DMMG and the PetscViewer");
478:   }

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

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

519:     Collective on DMMG

521:     Input Parameter:
522: +   dmmg - the context
523: .   has_cnst - is the constant vector in the null space
524: .   n - number of null vectors (excluding the possible constant vector)
525: -   func - a function that fills an array of vectors with the null vectors (must be orthonormal), may be PETSC_NULL

527:     Level: advanced

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

531: @*/
532: PetscErrorCode DMMGSetNullSpace(DMMG *dmmg,PetscTruth has_cnst,PetscInt n,PetscErrorCode (*func)(DMMG,Vec[]))
533: {
535:   PetscInt       i,j,nlevels = dmmg[0]->nlevels;
536:   Vec            *nulls = 0;
537:   MatNullSpace   nullsp;
538:   KSP            iksp;
539:   PC             pc,ipc;
540:   PetscTruth     ismg,isred;

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

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

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

591:     Collective on DMMG

593:     Input Parameter:
594: +   dmmg - the context
595: -   vec - dummy argument

597:     Level: intermediate

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

601: @*/
602: PetscErrorCode DMMGInitialGuessCurrent(DMMG dmmg,Vec vec)
603: {
605:   return(0);
606: }

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

613:     Collective on DMMG

615:     Input Parameter:
616: +   dmmg - the context
617: -   guess - the function

619:     Notes: For nonlinear problems, if this is not set, then the current value in the 
620:              solution vector (obtained with DMMGGetX()) is used. Thus is if you doing 'time
621:              stepping' it will use your current solution as the guess for the next timestep.
622:            If grid sequencing is used (via -dmmg_grid_sequence) then the "guess" function
623:              is used only on the coarsest grid.
624:            For linear problems, if this is not set, then 0 is used as an initial guess.
625:              If you would like the linear solver to also (like the nonlinear solver) use
626:              the current solution vector as the initial guess then use DMMGInitialGuessCurrent()
627:              as the function you pass in

629:     Level: intermediate


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

634: @*/
635: PetscErrorCode DMMGSetInitialGuess(DMMG *dmmg,PetscErrorCode (*guess)(DMMG,Vec))
636: {
637:   PetscInt i,nlevels = dmmg[0]->nlevels;

640:   for (i=0; i<nlevels; i++) {
641:     dmmg[i]->initialguess = guess;
642:   }
643:   return(0);
644: }