Actual source code: snesngmres.c

petsc-master 2019-06-15
Report Typos and Errors
  1:  #include <../src/snes/impls/ngmres/snesngmres.h>
  2:  #include <petscblaslapack.h>
  3:  #include <petscdm.h>

  5: const char *const SNESNGMRESRestartTypes[] = {"NONE","PERIODIC","DIFFERENCE","SNESNGMRESRestartType","SNES_NGMRES_RESTART_",0};
  6: const char *const SNESNGMRESSelectTypes[] = {"NONE","DIFFERENCE","LINESEARCH","SNESNGMRESSelectType","SNES_NGMRES_SELECT_",0};

  8: PetscErrorCode SNESReset_NGMRES(SNES snes)
  9: {
 10:   SNES_NGMRES    *ngmres = (SNES_NGMRES*) snes->data;

 14:   VecDestroyVecs(ngmres->msize,&ngmres->Fdot);
 15:   VecDestroyVecs(ngmres->msize,&ngmres->Xdot);
 16:   SNESLineSearchDestroy(&ngmres->additive_linesearch);
 17:   return(0);
 18: }

 20: PetscErrorCode SNESDestroy_NGMRES(SNES snes)
 21: {
 23:   SNES_NGMRES    *ngmres = (SNES_NGMRES*)snes->data;

 26:   SNESReset_NGMRES(snes);
 27:   PetscFree5(ngmres->h,ngmres->beta,ngmres->xi,ngmres->fnorms,ngmres->q);
 28:   PetscFree(ngmres->s);
 29:   PetscFree(ngmres->xnorms);
 30: #if defined(PETSC_USE_COMPLEX)
 31:   PetscFree(ngmres->rwork);
 32: #endif
 33:   PetscFree(ngmres->work);
 34:   PetscFree(snes->data);
 35:   return(0);
 36: }

 38: PetscErrorCode SNESSetUp_NGMRES(SNES snes)
 39: {
 40:   SNES_NGMRES    *ngmres = (SNES_NGMRES*) snes->data;
 41:   const char     *optionsprefix;
 42:   PetscInt       msize,hsize;
 44:   DM             dm;

 47:   if (snes->npc && snes->npcside== PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
 48:     SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE,"SNESNGMRES does not support left preconditioning with unpreconditioned function");
 49:   }
 50:   if (snes->npcside== PC_LEFT && snes->functype == SNES_FUNCTION_DEFAULT) snes->functype = SNES_FUNCTION_PRECONDITIONED;
 51:   SNESSetWorkVecs(snes,5);

 53:   if (!snes->vec_sol) {
 54:     SNESGetDM(snes,&dm);
 55:     DMCreateGlobalVector(dm,&snes->vec_sol);
 56:   }

 58:   if (!ngmres->Xdot) {VecDuplicateVecs(snes->vec_sol,ngmres->msize,&ngmres->Xdot);}
 59:   if (!ngmres->Fdot) {VecDuplicateVecs(snes->vec_sol,ngmres->msize,&ngmres->Fdot);}
 60:   if (!ngmres->setup_called) {
 61:     msize = ngmres->msize;          /* restart size */
 62:     hsize = msize * msize;

 64:     /* explicit least squares minimization solve */
 65:     PetscMalloc5(hsize,&ngmres->h, msize,&ngmres->beta, msize,&ngmres->xi, msize,&ngmres->fnorms, hsize,&ngmres->q);
 66:     PetscMalloc1(msize,&ngmres->xnorms);
 67:     ngmres->nrhs  = 1;
 68:     ngmres->lda   = msize;
 69:     ngmres->ldb   = msize;
 70:     PetscMalloc1(msize,&ngmres->s);
 71:     PetscMemzero(ngmres->h,   hsize*sizeof(PetscScalar));
 72:     PetscMemzero(ngmres->q,   hsize*sizeof(PetscScalar));
 73:     PetscMemzero(ngmres->xi,  msize*sizeof(PetscScalar));
 74:     PetscMemzero(ngmres->beta,msize*sizeof(PetscScalar));
 75:     ngmres->lwork = 12*msize;
 76: #if defined(PETSC_USE_COMPLEX)
 77:     PetscMalloc1(ngmres->lwork,&ngmres->rwork);
 78: #endif
 79:     PetscMalloc1(ngmres->lwork,&ngmres->work);
 80:   }

 82:   /* linesearch setup */
 83:   SNESGetOptionsPrefix(snes,&optionsprefix);

 85:   if (ngmres->select_type == SNES_NGMRES_SELECT_LINESEARCH) {
 86:     SNESLineSearchCreate(PetscObjectComm((PetscObject)snes),&ngmres->additive_linesearch);
 87:     SNESLineSearchSetSNES(ngmres->additive_linesearch,snes);
 88:     SNESLineSearchSetType(ngmres->additive_linesearch,SNESLINESEARCHL2);
 89:     SNESLineSearchAppendOptionsPrefix(ngmres->additive_linesearch,"additive_");
 90:     SNESLineSearchAppendOptionsPrefix(ngmres->additive_linesearch,optionsprefix);
 91:     SNESLineSearchSetFromOptions(ngmres->additive_linesearch);
 92:   }

 94:   ngmres->setup_called = PETSC_TRUE;
 95:   return(0);
 96: }

 98: PetscErrorCode SNESSetFromOptions_NGMRES(PetscOptionItems *PetscOptionsObject,SNES snes)
 99: {
100:   SNES_NGMRES    *ngmres = (SNES_NGMRES*) snes->data;
102:   PetscBool      debug = PETSC_FALSE;

105:   PetscOptionsHead(PetscOptionsObject,"SNES NGMRES options");
106:   PetscOptionsEnum("-snes_ngmres_select_type","Select type","SNESNGMRESSetSelectType",SNESNGMRESSelectTypes,
107:                           (PetscEnum)ngmres->select_type,(PetscEnum*)&ngmres->select_type,NULL);
108:   PetscOptionsEnum("-snes_ngmres_restart_type","Restart type","SNESNGMRESSetRestartType",SNESNGMRESRestartTypes,
109:                           (PetscEnum)ngmres->restart_type,(PetscEnum*)&ngmres->restart_type,NULL);
110:   PetscOptionsBool("-snes_ngmres_candidate", "Use candidate storage",              "SNES",ngmres->candidate,&ngmres->candidate,NULL);
111:   PetscOptionsBool("-snes_ngmres_approxfunc","Linearly approximate the function", "SNES",ngmres->approxfunc,&ngmres->approxfunc,NULL);
112:   PetscOptionsInt("-snes_ngmres_m",          "Number of directions",               "SNES",ngmres->msize,&ngmres->msize,NULL);
113:   PetscOptionsInt("-snes_ngmres_restart",    "Iterations before forced restart",   "SNES",ngmres->restart_periodic,&ngmres->restart_periodic,NULL);
114:   PetscOptionsInt("-snes_ngmres_restart_it", "Tolerance iterations before restart","SNES",ngmres->restart_it,&ngmres->restart_it,NULL);
115:   PetscOptionsBool("-snes_ngmres_monitor",   "Monitor actions of NGMRES",          "SNES",ngmres->monitor ? PETSC_TRUE : PETSC_FALSE,&debug,NULL);
116:   if (debug) {
117:     ngmres->monitor = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes));
118:   }
119:   PetscOptionsReal("-snes_ngmres_gammaA",    "Residual selection constant",   "SNES",ngmres->gammaA,&ngmres->gammaA,NULL);
120:   PetscOptionsReal("-snes_ngmres_gammaC",    "Residual restart constant",     "SNES",ngmres->gammaC,&ngmres->gammaC,NULL);
121:   PetscOptionsReal("-snes_ngmres_epsilonB",  "Difference selection constant", "SNES",ngmres->epsilonB,&ngmres->epsilonB,NULL);
122:   PetscOptionsReal("-snes_ngmres_deltaB",    "Difference residual selection constant", "SNES",ngmres->deltaB,&ngmres->deltaB,NULL);
123:   PetscOptionsBool("-snes_ngmres_single_reduction", "Aggregate reductions",  "SNES",ngmres->singlereduction,&ngmres->singlereduction,NULL);
124:   PetscOptionsBool("-snes_ngmres_restart_fm_rise", "Restart on F_M residual rise",  "SNESNGMRESSetRestartFmRise",ngmres->restart_fm_rise,&ngmres->restart_fm_rise,NULL);
125:   PetscOptionsTail();
126:   if ((ngmres->gammaA > ngmres->gammaC) && (ngmres->gammaC > 2.)) ngmres->gammaC = ngmres->gammaA;
127:   return(0);
128: }

130: PetscErrorCode SNESView_NGMRES(SNES snes,PetscViewer viewer)
131: {
132:   SNES_NGMRES    *ngmres = (SNES_NGMRES*) snes->data;
133:   PetscBool      iascii;

137:   PetscObjectTypeCompare((PetscObject) viewer,PETSCVIEWERASCII,&iascii);
138:   if (iascii) {
139:     PetscViewerASCIIPrintf(viewer,"  Number of stored past updates: %d\n", ngmres->msize);
140:     PetscViewerASCIIPrintf(viewer,"  Residual selection: gammaA=%1.0e, gammaC=%1.0e\n",ngmres->gammaA,ngmres->gammaC);
141:     PetscViewerASCIIPrintf(viewer,"  Difference restart: epsilonB=%1.0e, deltaB=%1.0e\n",ngmres->epsilonB,ngmres->deltaB);
142:     PetscViewerASCIIPrintf(viewer,"  Restart on F_M residual increase: %s\n",ngmres->restart_fm_rise?"TRUE":"FALSE");
143:   }
144:   return(0);
145: }

147: PetscErrorCode SNESSolve_NGMRES(SNES snes)
148: {
149:   SNES_NGMRES          *ngmres = (SNES_NGMRES*) snes->data;
150:   /* present solution, residual, and preconditioned residual */
151:   Vec                  X,F,B,D,Y;

153:   /* candidate linear combination answers */
154:   Vec                  XA,FA,XM,FM;

156:   /* coefficients and RHS to the minimization problem */
157:   PetscReal            fnorm,fMnorm,fAnorm;
158:   PetscReal            xnorm,xMnorm,xAnorm;
159:   PetscReal            ynorm,yMnorm,yAnorm;
160:   PetscInt             k,k_restart,l,ivec,restart_count = 0;

162:   /* solution selection data */
163:   PetscBool            selectRestart;
164:   /*
165:       These two variables are initialized to prevent compilers/analyzers from producing false warnings about these variables being passed
166:       to SNESNGMRESSelect_Private() without being set when SNES_NGMRES_RESTART_DIFFERENCE, the values are not used in the subroutines in that case
167:       so the code is correct as written.
168:   */
169:   PetscReal            dnorm = 0.0,dminnorm = 0.0;
170:   PetscReal            fminnorm;

172:   SNESConvergedReason  reason;
173:   SNESLineSearchReason lssucceed;
174:   PetscErrorCode       ierr;

177:   if (snes->xl || snes->xu || snes->ops->computevariablebounds) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name);

179:   PetscCitationsRegister(SNESCitation,&SNEScite);
180:   /* variable initialization */
181:   snes->reason = SNES_CONVERGED_ITERATING;
182:   X            = snes->vec_sol;
183:   F            = snes->vec_func;
184:   B            = snes->vec_rhs;
185:   XA           = snes->vec_sol_update;
186:   FA           = snes->work[0];
187:   D            = snes->work[1];

189:   /* work for the line search */
190:   Y  = snes->work[2];
191:   XM = snes->work[3];
192:   FM = snes->work[4];

194:   PetscObjectSAWsTakeAccess((PetscObject)snes);
195:   snes->iter = 0;
196:   snes->norm = 0.;
197:   PetscObjectSAWsGrantAccess((PetscObject)snes);

199:   /* initialization */

201:   if (snes->npc && snes->npcside== PC_LEFT) {
202:     SNESApplyNPC(snes,X,NULL,F);
203:     SNESGetConvergedReason(snes->npc,&reason);
204:     if (reason < 0  && reason != SNES_DIVERGED_MAX_IT) {
205:       snes->reason = SNES_DIVERGED_INNER;
206:       return(0);
207:     }
208:     VecNorm(F,NORM_2,&fnorm);
209:   } else {
210:     if (!snes->vec_func_init_set) {
211:       SNESComputeFunction(snes,X,F);
212:     } else snes->vec_func_init_set = PETSC_FALSE;

214:     VecNorm(F,NORM_2,&fnorm);
215:     SNESCheckFunctionNorm(snes,fnorm);
216:   }
217:   fminnorm = fnorm;

219:   PetscObjectSAWsTakeAccess((PetscObject)snes);
220:   snes->norm = fnorm;
221:   PetscObjectSAWsGrantAccess((PetscObject)snes);
222:   SNESLogConvergenceHistory(snes,fnorm,0);
223:   SNESMonitor(snes,0,fnorm);
224:   (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);
225:   if (snes->reason) return(0);
226:   SNESNGMRESUpdateSubspace_Private(snes,0,0,F,fnorm,X);

228:   k_restart = 1;
229:   l         = 1;
230:   ivec      = 0;
231:   for (k=1; k < snes->max_its+1; k++) {
232:     /* Computation of x^M */
233:     if (snes->npc && snes->npcside== PC_RIGHT) {
234:       VecCopy(X,XM);
235:       SNESSetInitialFunction(snes->npc,F);

237:       PetscLogEventBegin(SNES_NPCSolve,snes->npc,XM,B,0);
238:       SNESSolve(snes->npc,B,XM);
239:       PetscLogEventEnd(SNES_NPCSolve,snes->npc,XM,B,0);

241:       SNESGetConvergedReason(snes->npc,&reason);
242:       if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
243:         snes->reason = SNES_DIVERGED_INNER;
244:         return(0);
245:       }
246:       SNESGetNPCFunction(snes,FM,&fMnorm);
247:     } else {
248:       /* no preconditioner -- just take gradient descent with line search */
249:       VecCopy(F,Y);
250:       VecCopy(F,FM);
251:       VecCopy(X,XM);

253:       fMnorm = fnorm;

255:       SNESLineSearchApply(snes->linesearch,XM,FM,&fMnorm,Y);
256:       SNESLineSearchGetReason(snes->linesearch,&lssucceed);
257:       if (lssucceed) {
258:         if (++snes->numFailures >= snes->maxFailures) {
259:           snes->reason = SNES_DIVERGED_LINE_SEARCH;
260:           return(0);
261:         }
262:       }
263:     }

265:     SNESNGMRESFormCombinedSolution_Private(snes,ivec,l,XM,FM,fMnorm,X,XA,FA);
266:     /* r = F(x) */
267:     if (fminnorm > fMnorm) fminnorm = fMnorm;  /* the minimum norm is now of F^M */

269:     /* differences for selection and restart */
270:     if (ngmres->restart_type == SNES_NGMRES_RESTART_DIFFERENCE || ngmres->select_type == SNES_NGMRES_SELECT_DIFFERENCE) {
271:       SNESNGMRESNorms_Private(snes,l,X,F,XM,FM,XA,FA,D,&dnorm,&dminnorm,&xMnorm,NULL,&yMnorm,&xAnorm,&fAnorm,&yAnorm);
272:     } else {
273:       SNESNGMRESNorms_Private(snes,l,X,F,XM,FM,XA,FA,D,NULL,NULL,&xMnorm,NULL,&yMnorm,&xAnorm,&fAnorm,&yAnorm);
274:     }
275:     SNESCheckFunctionNorm(snes,fnorm);

277:     /* combination (additive) or selection (multiplicative) of the N-GMRES solution */
278:     SNESNGMRESSelect_Private(snes,k_restart,XM,FM,xMnorm,fMnorm,yMnorm,XA,FA,xAnorm,fAnorm,yAnorm,dnorm,fminnorm,dminnorm,X,F,Y,&xnorm,&fnorm,&ynorm);
279:     selectRestart = PETSC_FALSE;

281:     if (ngmres->restart_type == SNES_NGMRES_RESTART_DIFFERENCE) {
282:       SNESNGMRESSelectRestart_Private(snes,l,fMnorm,fAnorm,dnorm,fminnorm,dminnorm,&selectRestart);

284:       /* if the restart conditions persist for more than restart_it iterations, restart. */
285:       if (selectRestart) restart_count++;
286:       else restart_count = 0;
287:     } else if (ngmres->restart_type == SNES_NGMRES_RESTART_PERIODIC) {
288:       if (k_restart > ngmres->restart_periodic) {
289:         if (ngmres->monitor) PetscViewerASCIIPrintf(ngmres->monitor,"periodic restart after %D iterations\n",k_restart);
290:         restart_count = ngmres->restart_it;
291:       }
292:     }

294:     ivec = k_restart % ngmres->msize; /* replace the last used part of the subspace */

296:     /* restart after restart conditions have persisted for a fixed number of iterations */
297:     if (restart_count >= ngmres->restart_it) {
298:       if (ngmres->monitor) {
299:         PetscViewerASCIIPrintf(ngmres->monitor,"Restarted at iteration %d\n",k_restart);
300:       }
301:       restart_count = 0;
302:       k_restart     = 1;
303:       l             = 1;
304:       ivec          = 0;
305:       /* q_{00} = nu */
306:       SNESNGMRESUpdateSubspace_Private(snes,0,0,FM,fMnorm,XM);
307:     } else {
308:       /* select the current size of the subspace */
309:       if (l < ngmres->msize) l++;
310:       k_restart++;
311:       /* place the current entry in the list of previous entries */
312:       if (ngmres->candidate) {
313:         if (fminnorm > fMnorm) fminnorm = fMnorm;
314:         SNESNGMRESUpdateSubspace_Private(snes,ivec,l,FM,fMnorm,XM);
315:       } else {
316:         if (fminnorm > fnorm) fminnorm = fnorm;
317:         SNESNGMRESUpdateSubspace_Private(snes,ivec,l,F,fnorm,X);
318:       }
319:     }

321:     PetscObjectSAWsTakeAccess((PetscObject)snes);
322:     snes->iter = k;
323:     snes->norm = fnorm;
324:     PetscObjectSAWsGrantAccess((PetscObject)snes);
325:     SNESLogConvergenceHistory(snes,snes->norm,snes->iter);
326:     SNESMonitor(snes,snes->iter,snes->norm);
327:     (*snes->ops->converged)(snes,snes->iter,0,0,fnorm,&snes->reason,snes->cnvP);
328:     if (snes->reason) return(0);
329:   }
330:   snes->reason = SNES_DIVERGED_MAX_IT;
331:   return(0);
332: }

334: /*@
335:  SNESNGMRESSetRestartFmRise - Increase the restart count if the step x_M increases the residual F_M

337:   Input Parameters:
338:   +  snes - the SNES context.
339:   -  flg  - boolean value deciding whether to use the option or not

341:   Options Database:
342:   + -snes_ngmres_restart_fm_rise - Increase the restart count if the step x_M increases the residual F_M

344:   Level: intermediate

346:   Notes:
347:   If the proposed step x_M increases the residual F_M, it might be trying to get out of a stagnation area.
348:   To help the solver do that, reset the Krylov subspace whenever F_M increases.

350:   This option must be used with SNES_NGMRES_RESTART_DIFFERENCE

352:   The default is FALSE.
353:   .seealso: SNES_NGMRES_RESTART_DIFFERENCE
354:   @*/
355: PetscErrorCode SNESNGMRESSetRestartFmRise(SNES snes,PetscBool flg)
356: {
357:     PetscErrorCode (*f)(SNES,PetscBool);

361:     PetscObjectQueryFunction((PetscObject)snes,"SNESNGMRESSetRestartFmRise_C",&f);
362:     if (f) {(f)(snes,flg);}
363:     return(0);
364: }

366: PetscErrorCode SNESNGMRESSetRestartFmRise_NGMRES(SNES snes,PetscBool flg)
367: {
368:   SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data;

371:   ngmres->restart_fm_rise = flg;
372:   return(0);
373: }

375: PetscErrorCode SNESNGMRESGetRestartFmRise(SNES snes,PetscBool *flg)
376: {
377:     PetscErrorCode (*f)(SNES,PetscBool*);

381:     PetscObjectQueryFunction((PetscObject)snes,"SNESNGMRESGetRestartFmRise_C",&f);
382:     if (f) {(f)(snes,flg);}
383:     return(0);
384: }

386: PetscErrorCode SNESNGMRESGetRestartFmRise_NGMRES(SNES snes,PetscBool *flg)
387: {
388:   SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data;

391:   *flg = ngmres->restart_fm_rise;
392:   return(0);
393: }


396: /*@
397:     SNESNGMRESSetRestartType - Sets the restart type for SNESNGMRES.

399:     Logically Collective on SNES

401:     Input Parameters:
402: +   snes - the iterative context
403: -   rtype - restart type

405:     Options Database:
406: +   -snes_ngmres_restart_type<difference,periodic,none> - set the restart type
407: -   -snes_ngmres_restart[30] - sets the number of iterations before restart for periodic

409:     Level: intermediate

411:     SNESNGMRESRestartTypes:
412: +   SNES_NGMRES_RESTART_NONE - never restart
413: .   SNES_NGMRES_RESTART_DIFFERENCE - restart based upon difference criteria
414: -   SNES_NGMRES_RESTART_PERIODIC - restart after a fixed number of iterations

416:     Notes:
417:     The default line search used is the L2 line search and it requires two additional function evaluations.

419: @*/
420: PetscErrorCode SNESNGMRESSetRestartType(SNES snes,SNESNGMRESRestartType rtype)
421: {

426:   PetscTryMethod(snes,"SNESNGMRESSetRestartType_C",(SNES,SNESNGMRESRestartType),(snes,rtype));
427:   return(0);
428: }

430: /*@
431:     SNESNGMRESSetSelectType - Sets the selection type for SNESNGMRES.  This determines how the candidate solution and
432:     combined solution are used to create the next iterate.

434:     Logically Collective on SNES

436:     Input Parameters:
437: +   snes - the iterative context
438: -   stype - selection type

440:     Options Database:
441: .   -snes_ngmres_select_type<difference,none,linesearch>

443:     Level: intermediate

445:     SNESNGMRESSelectTypes:
446: +   SNES_NGMRES_SELECT_NONE - choose the combined solution all the time
447: .   SNES_NGMRES_SELECT_DIFFERENCE - choose based upon the selection criteria
448: -   SNES_NGMRES_SELECT_LINESEARCH - choose based upon line search combination

450:     Notes:
451:     The default line search used is the L2 line search and it requires two additional function evaluations.

453: @*/
454: PetscErrorCode SNESNGMRESSetSelectType(SNES snes,SNESNGMRESSelectType stype)
455: {

460:   PetscTryMethod(snes,"SNESNGMRESSetSelectType_C",(SNES,SNESNGMRESSelectType),(snes,stype));
461:   return(0);
462: }

464: PetscErrorCode SNESNGMRESSetSelectType_NGMRES(SNES snes,SNESNGMRESSelectType stype)
465: {
466:   SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data;

469:   ngmres->select_type = stype;
470:   return(0);
471: }

473: PetscErrorCode SNESNGMRESSetRestartType_NGMRES(SNES snes,SNESNGMRESRestartType rtype)
474: {
475:   SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data;

478:   ngmres->restart_type = rtype;
479:   return(0);
480: }

482: /*MC
483:   SNESNGMRES - The Nonlinear Generalized Minimum Residual method.

485:    Level: beginner

487:    Options Database:
488: +  -snes_ngmres_select_type<difference,none,linesearch> - choose the select between candidate and combined solution
489: .  -snes_ngmres_restart_type<difference,none,periodic> - choose the restart conditions
490: .  -snes_ngmres_candidate        - Use NGMRES variant which combines candidate solutions instead of actual solutions
491: .  -snes_ngmres_m                - Number of stored previous solutions and residuals
492: .  -snes_ngmres_restart_it       - Number of iterations the restart conditions hold before restart
493: .  -snes_ngmres_gammaA           - Residual tolerance for solution select between the candidate and combination
494: .  -snes_ngmres_gammaC           - Residual tolerance for restart
495: .  -snes_ngmres_epsilonB         - Difference tolerance between subsequent solutions triggering restart
496: .  -snes_ngmres_deltaB           - Difference tolerance between residuals triggering restart
497: .  -snes_ngmres_restart_fm_rise  - Restart on residual rise from x_M step
498: .  -snes_ngmres_monitor          - Prints relevant information about the ngmres iteration
499: .  -snes_linesearch_type <basic,l2,cp> - Line search type used for the default smoother
500: -  -additive_snes_linesearch_type - linesearch type used to select between the candidate and combined solution with additive select type

502:    Notes:

504:    The N-GMRES method combines m previous solutions into a minimum-residual solution by solving a small linearized
505:    optimization problem at each iteration.

507:    Very similar to the SNESANDERSON algorithm.

509:    References:
510: +  1. - C. W. Oosterlee and T. Washio, "Krylov Subspace Acceleration of Nonlinear Multigrid with Application to Recirculating Flows", 
511:    SIAM Journal on Scientific Computing, 21(5), 2000.
512: -  2. - Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu, "Composing Scalable Nonlinear Algebraic Solvers", 
513:    SIAM Review, 57(4), 2015


516: .seealso: SNESCreate(), SNES, SNESSetType(), SNESType (for list of available types)
517: M*/

519: PETSC_EXTERN PetscErrorCode SNESCreate_NGMRES(SNES snes)
520: {
521:   SNES_NGMRES    *ngmres;
523:   SNESLineSearch linesearch;

526:   snes->ops->destroy        = SNESDestroy_NGMRES;
527:   snes->ops->setup          = SNESSetUp_NGMRES;
528:   snes->ops->setfromoptions = SNESSetFromOptions_NGMRES;
529:   snes->ops->view           = SNESView_NGMRES;
530:   snes->ops->solve          = SNESSolve_NGMRES;
531:   snes->ops->reset          = SNESReset_NGMRES;

533:   snes->usesnpc  = PETSC_TRUE;
534:   snes->usesksp  = PETSC_FALSE;
535:   snes->npcside  = PC_RIGHT;

537:   snes->alwayscomputesfinalresidual = PETSC_TRUE;

539:   PetscNewLog(snes,&ngmres);
540:   snes->data    = (void*) ngmres;
541:   ngmres->msize = 30;

543:   if (!snes->tolerancesset) {
544:     snes->max_funcs = 30000;
545:     snes->max_its   = 10000;
546:   }

548:   ngmres->candidate = PETSC_FALSE;

550:  SNESGetLineSearch(snes,&linesearch);
551:  SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC);

553:   ngmres->additive_linesearch = NULL;
554:   ngmres->approxfunc          = PETSC_FALSE;
555:   ngmres->restart_it          = 2;
556:   ngmres->restart_periodic    = 30;
557:   ngmres->gammaA              = 2.0;
558:   ngmres->gammaC              = 2.0;
559:   ngmres->deltaB              = 0.9;
560:   ngmres->epsilonB            = 0.1;
561:   ngmres->restart_fm_rise     = PETSC_FALSE;

563:   ngmres->restart_type = SNES_NGMRES_RESTART_DIFFERENCE;
564:   ngmres->select_type  = SNES_NGMRES_SELECT_DIFFERENCE;

566:   PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetSelectType_C",SNESNGMRESSetSelectType_NGMRES);
567:   PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetRestartType_C",SNESNGMRESSetRestartType_NGMRES);
568:   PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetRestartFmRise_C",SNESNGMRESSetRestartFmRise_NGMRES);
569:   PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESGetRestartFmRise_C",SNESNGMRESGetRestartFmRise_NGMRES);
570:   return(0);
571: }