xref: /petsc/src/snes/impls/ngmres/snesngmres.c (revision 117ef88edefbfc12e7c19efe87a19a2e1b0acd4f)
1 #include <../src/snes/impls/ngmres/snesngmres.h> /*I "petscsnes.h" I*/
2 #include <petscblaslapack.h>
3 #include <petscdm.h>
4 
5 const char *const SNESNGMRESRestartTypes[] = {"NONE","PERIODIC","DIFFERENCE","SNESNGMRESRestartType","SNES_NGMRES_RESTART_",NULL};
6 const char *const SNESNGMRESSelectTypes[] = {"NONE","DIFFERENCE","LINESEARCH","SNESNGMRESSelectType","SNES_NGMRES_SELECT_",NULL};
7 
8 PetscErrorCode SNESReset_NGMRES(SNES snes)
9 {
10   SNES_NGMRES    *ngmres = (SNES_NGMRES*) snes->data;
11   PetscErrorCode ierr;
12 
13   PetscFunctionBegin;
14   ierr = VecDestroyVecs(ngmres->msize,&ngmres->Fdot);CHKERRQ(ierr);
15   ierr = VecDestroyVecs(ngmres->msize,&ngmres->Xdot);CHKERRQ(ierr);
16   ierr = SNESLineSearchDestroy(&ngmres->additive_linesearch);CHKERRQ(ierr);
17   PetscFunctionReturn(0);
18 }
19 
20 PetscErrorCode SNESDestroy_NGMRES(SNES snes)
21 {
22   PetscErrorCode ierr;
23   SNES_NGMRES    *ngmres = (SNES_NGMRES*)snes->data;
24 
25   PetscFunctionBegin;
26   ierr = SNESReset_NGMRES(snes);CHKERRQ(ierr);
27   ierr = PetscFree4(ngmres->h,ngmres->beta,ngmres->xi,ngmres->q);CHKERRQ(ierr);
28   ierr = PetscFree3(ngmres->xnorms,ngmres->fnorms,ngmres->s);CHKERRQ(ierr);
29 #if defined(PETSC_USE_COMPLEX)
30   ierr = PetscFree(ngmres->rwork);CHKERRQ(ierr);
31 #endif
32   ierr = PetscFree(ngmres->work);CHKERRQ(ierr);
33   ierr = PetscFree(snes->data);CHKERRQ(ierr);
34   PetscFunctionReturn(0);
35 }
36 
37 PetscErrorCode SNESSetUp_NGMRES(SNES snes)
38 {
39   SNES_NGMRES    *ngmres = (SNES_NGMRES*) snes->data;
40   const char     *optionsprefix;
41   PetscInt       msize,hsize;
42   PetscErrorCode ierr;
43   DM             dm;
44 
45   PetscFunctionBegin;
46   if (snes->npc && snes->npcside== PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
47     SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE,"SNESNGMRES does not support left preconditioning with unpreconditioned function");
48   }
49   if (snes->npcside== PC_LEFT && snes->functype == SNES_FUNCTION_DEFAULT) snes->functype = SNES_FUNCTION_PRECONDITIONED;
50   ierr = SNESSetWorkVecs(snes,5);CHKERRQ(ierr);
51 
52   if (!snes->vec_sol) {
53     ierr             = SNESGetDM(snes,&dm);CHKERRQ(ierr);
54     ierr             = DMCreateGlobalVector(dm,&snes->vec_sol);CHKERRQ(ierr);
55   }
56 
57   if (!ngmres->Xdot) {ierr = VecDuplicateVecs(snes->vec_sol,ngmres->msize,&ngmres->Xdot);CHKERRQ(ierr);}
58   if (!ngmres->Fdot) {ierr = VecDuplicateVecs(snes->vec_sol,ngmres->msize,&ngmres->Fdot);CHKERRQ(ierr);}
59   if (!ngmres->setup_called) {
60     msize = ngmres->msize;          /* restart size */
61     hsize = msize * msize;
62 
63     /* explicit least squares minimization solve */
64     ierr = PetscCalloc4(hsize,&ngmres->h, msize,&ngmres->beta, msize,&ngmres->xi, hsize,&ngmres->q);CHKERRQ(ierr);
65     ierr = PetscMalloc3(msize,&ngmres->xnorms,msize,&ngmres->fnorms,msize,&ngmres->s);CHKERRQ(ierr);
66     ngmres->nrhs  = 1;
67     ngmres->lda   = msize;
68     ngmres->ldb   = msize;
69     ngmres->lwork = 12*msize;
70 #if defined(PETSC_USE_COMPLEX)
71     ierr = PetscMalloc1(ngmres->lwork,&ngmres->rwork);CHKERRQ(ierr);
72 #endif
73     ierr = PetscMalloc1(ngmres->lwork,&ngmres->work);CHKERRQ(ierr);
74   }
75 
76   /* linesearch setup */
77   ierr = SNESGetOptionsPrefix(snes,&optionsprefix);CHKERRQ(ierr);
78 
79   if (ngmres->select_type == SNES_NGMRES_SELECT_LINESEARCH) {
80     ierr = SNESLineSearchCreate(PetscObjectComm((PetscObject)snes),&ngmres->additive_linesearch);CHKERRQ(ierr);
81     ierr = SNESLineSearchSetSNES(ngmres->additive_linesearch,snes);CHKERRQ(ierr);
82     if (!((PetscObject)ngmres->additive_linesearch)->type_name) {
83       ierr = SNESLineSearchSetType(ngmres->additive_linesearch,SNESLINESEARCHL2);CHKERRQ(ierr);
84     }
85     ierr = SNESLineSearchAppendOptionsPrefix(ngmres->additive_linesearch,"additive_");CHKERRQ(ierr);
86     ierr = SNESLineSearchAppendOptionsPrefix(ngmres->additive_linesearch,optionsprefix);CHKERRQ(ierr);
87     ierr = SNESLineSearchSetFromOptions(ngmres->additive_linesearch);CHKERRQ(ierr);
88   }
89 
90   ngmres->setup_called = PETSC_TRUE;
91   PetscFunctionReturn(0);
92 }
93 
94 PetscErrorCode SNESSetFromOptions_NGMRES(PetscOptionItems *PetscOptionsObject,SNES snes)
95 {
96   SNES_NGMRES    *ngmres = (SNES_NGMRES*) snes->data;
97   PetscErrorCode ierr;
98   PetscBool      debug = PETSC_FALSE;
99 
100   PetscFunctionBegin;
101   ierr = PetscOptionsHead(PetscOptionsObject,"SNES NGMRES options");CHKERRQ(ierr);
102   ierr = PetscOptionsEnum("-snes_ngmres_select_type","Select type","SNESNGMRESSetSelectType",SNESNGMRESSelectTypes,
103                           (PetscEnum)ngmres->select_type,(PetscEnum*)&ngmres->select_type,NULL);CHKERRQ(ierr);
104   ierr = PetscOptionsEnum("-snes_ngmres_restart_type","Restart type","SNESNGMRESSetRestartType",SNESNGMRESRestartTypes,
105                           (PetscEnum)ngmres->restart_type,(PetscEnum*)&ngmres->restart_type,NULL);CHKERRQ(ierr);
106   ierr = PetscOptionsBool("-snes_ngmres_candidate", "Use candidate storage",              "SNES",ngmres->candidate,&ngmres->candidate,NULL);CHKERRQ(ierr);
107   ierr = PetscOptionsBool("-snes_ngmres_approxfunc","Linearly approximate the function", "SNES",ngmres->approxfunc,&ngmres->approxfunc,NULL);CHKERRQ(ierr);
108   ierr = PetscOptionsInt("-snes_ngmres_m",          "Number of directions",               "SNES",ngmres->msize,&ngmres->msize,NULL);CHKERRQ(ierr);
109   ierr = PetscOptionsInt("-snes_ngmres_restart",    "Iterations before forced restart",   "SNES",ngmres->restart_periodic,&ngmres->restart_periodic,NULL);CHKERRQ(ierr);
110   ierr = PetscOptionsInt("-snes_ngmres_restart_it", "Tolerance iterations before restart","SNES",ngmres->restart_it,&ngmres->restart_it,NULL);CHKERRQ(ierr);
111   ierr = PetscOptionsBool("-snes_ngmres_monitor",   "Monitor actions of NGMRES",          "SNES",ngmres->monitor ? PETSC_TRUE : PETSC_FALSE,&debug,NULL);CHKERRQ(ierr);
112   if (debug) {
113     ngmres->monitor = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes));CHKERRQ(ierr);
114   }
115   ierr = PetscOptionsReal("-snes_ngmres_gammaA",    "Residual selection constant",   "SNES",ngmres->gammaA,&ngmres->gammaA,NULL);CHKERRQ(ierr);
116   ierr = PetscOptionsReal("-snes_ngmres_gammaC",    "Residual restart constant",     "SNES",ngmres->gammaC,&ngmres->gammaC,NULL);CHKERRQ(ierr);
117   ierr = PetscOptionsReal("-snes_ngmres_epsilonB",  "Difference selection constant", "SNES",ngmres->epsilonB,&ngmres->epsilonB,NULL);CHKERRQ(ierr);
118   ierr = PetscOptionsReal("-snes_ngmres_deltaB",    "Difference residual selection constant", "SNES",ngmres->deltaB,&ngmres->deltaB,NULL);CHKERRQ(ierr);
119   ierr = PetscOptionsBool("-snes_ngmres_single_reduction", "Aggregate reductions",  "SNES",ngmres->singlereduction,&ngmres->singlereduction,NULL);CHKERRQ(ierr);
120   ierr = PetscOptionsBool("-snes_ngmres_restart_fm_rise", "Restart on F_M residual rise",  "SNESNGMRESSetRestartFmRise",ngmres->restart_fm_rise,&ngmres->restart_fm_rise,NULL);CHKERRQ(ierr);
121   ierr = PetscOptionsTail();CHKERRQ(ierr);
122   if ((ngmres->gammaA > ngmres->gammaC) && (ngmres->gammaC > 2.)) ngmres->gammaC = ngmres->gammaA;
123   PetscFunctionReturn(0);
124 }
125 
126 PetscErrorCode SNESView_NGMRES(SNES snes,PetscViewer viewer)
127 {
128   SNES_NGMRES    *ngmres = (SNES_NGMRES*) snes->data;
129   PetscBool      iascii;
130   PetscErrorCode ierr;
131 
132   PetscFunctionBegin;
133   ierr = PetscObjectTypeCompare((PetscObject) viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr);
134   if (iascii) {
135     ierr = PetscViewerASCIIPrintf(viewer,"  Number of stored past updates: %d\n", ngmres->msize);CHKERRQ(ierr);
136     ierr = PetscViewerASCIIPrintf(viewer,"  Residual selection: gammaA=%1.0e, gammaC=%1.0e\n",ngmres->gammaA,ngmres->gammaC);CHKERRQ(ierr);
137     ierr = PetscViewerASCIIPrintf(viewer,"  Difference restart: epsilonB=%1.0e, deltaB=%1.0e\n",ngmres->epsilonB,ngmres->deltaB);CHKERRQ(ierr);
138     ierr = PetscViewerASCIIPrintf(viewer,"  Restart on F_M residual increase: %s\n",ngmres->restart_fm_rise?"TRUE":"FALSE");CHKERRQ(ierr);
139   }
140   PetscFunctionReturn(0);
141 }
142 
143 PetscErrorCode SNESSolve_NGMRES(SNES snes)
144 {
145   SNES_NGMRES          *ngmres = (SNES_NGMRES*) snes->data;
146   /* present solution, residual, and preconditioned residual */
147   Vec                  X,F,B,D,Y;
148 
149   /* candidate linear combination answers */
150   Vec                  XA,FA,XM,FM;
151 
152   /* coefficients and RHS to the minimization problem */
153   PetscReal            fnorm,fMnorm,fAnorm;
154   PetscReal            xnorm,xMnorm,xAnorm;
155   PetscReal            ynorm,yMnorm,yAnorm;
156   PetscInt             k,k_restart,l,ivec,restart_count = 0;
157 
158   /* solution selection data */
159   PetscBool            selectRestart;
160   /*
161       These two variables are initialized to prevent compilers/analyzers from producing false warnings about these variables being passed
162       to SNESNGMRESSelect_Private() without being set when SNES_NGMRES_RESTART_DIFFERENCE, the values are not used in the subroutines in that case
163       so the code is correct as written.
164   */
165   PetscReal            dnorm = 0.0,dminnorm = 0.0;
166   PetscReal            fminnorm;
167 
168   SNESConvergedReason  reason;
169   SNESLineSearchReason lssucceed;
170   PetscErrorCode       ierr;
171 
172   PetscFunctionBegin;
173   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);
174 
175   ierr = PetscCitationsRegister(SNESCitation,&SNEScite);CHKERRQ(ierr);
176   /* variable initialization */
177   snes->reason = SNES_CONVERGED_ITERATING;
178   X            = snes->vec_sol;
179   F            = snes->vec_func;
180   B            = snes->vec_rhs;
181   XA           = snes->vec_sol_update;
182   FA           = snes->work[0];
183   D            = snes->work[1];
184 
185   /* work for the line search */
186   Y  = snes->work[2];
187   XM = snes->work[3];
188   FM = snes->work[4];
189 
190   ierr       = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr);
191   snes->iter = 0;
192   snes->norm = 0.;
193   ierr       = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr);
194 
195   /* initialization */
196 
197   if (snes->npc && snes->npcside== PC_LEFT) {
198     ierr = SNESApplyNPC(snes,X,NULL,F);CHKERRQ(ierr);
199     ierr = SNESGetConvergedReason(snes->npc,&reason);CHKERRQ(ierr);
200     if (reason < 0  && reason != SNES_DIVERGED_MAX_IT) {
201       snes->reason = SNES_DIVERGED_INNER;
202       PetscFunctionReturn(0);
203     }
204     ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr);
205   } else {
206     if (!snes->vec_func_init_set) {
207       ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr);
208     } else snes->vec_func_init_set = PETSC_FALSE;
209 
210     ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr);
211     SNESCheckFunctionNorm(snes,fnorm);
212   }
213   fminnorm = fnorm;
214 
215   ierr       = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr);
216   snes->norm = fnorm;
217   ierr       = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr);
218   ierr       = SNESLogConvergenceHistory(snes,fnorm,0);CHKERRQ(ierr);
219   ierr       = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr);
220   ierr       = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr);
221   if (snes->reason) PetscFunctionReturn(0);
222   SNESNGMRESUpdateSubspace_Private(snes,0,0,F,fnorm,X);
223 
224   k_restart = 1;
225   l         = 1;
226   ivec      = 0;
227   for (k=1; k < snes->max_its+1; k++) {
228     /* Computation of x^M */
229     if (snes->npc && snes->npcside== PC_RIGHT) {
230       ierr = VecCopy(X,XM);CHKERRQ(ierr);
231       ierr = SNESSetInitialFunction(snes->npc,F);CHKERRQ(ierr);
232 
233       ierr = PetscLogEventBegin(SNES_NPCSolve,snes->npc,XM,B,0);CHKERRQ(ierr);
234       ierr = SNESSolve(snes->npc,B,XM);CHKERRQ(ierr);
235       ierr = PetscLogEventEnd(SNES_NPCSolve,snes->npc,XM,B,0);CHKERRQ(ierr);
236 
237       ierr = SNESGetConvergedReason(snes->npc,&reason);CHKERRQ(ierr);
238       if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
239         snes->reason = SNES_DIVERGED_INNER;
240         PetscFunctionReturn(0);
241       }
242       ierr = SNESGetNPCFunction(snes,FM,&fMnorm);CHKERRQ(ierr);
243     } else {
244       /* no preconditioner -- just take gradient descent with line search */
245       ierr = VecCopy(F,Y);CHKERRQ(ierr);
246       ierr = VecCopy(F,FM);CHKERRQ(ierr);
247       ierr = VecCopy(X,XM);CHKERRQ(ierr);
248 
249       fMnorm = fnorm;
250 
251       ierr = SNESLineSearchApply(snes->linesearch,XM,FM,&fMnorm,Y);CHKERRQ(ierr);
252       ierr = SNESLineSearchGetReason(snes->linesearch,&lssucceed);CHKERRQ(ierr);
253       if (lssucceed) {
254         if (++snes->numFailures >= snes->maxFailures) {
255           snes->reason = SNES_DIVERGED_LINE_SEARCH;
256           PetscFunctionReturn(0);
257         }
258       }
259     }
260 
261     ierr = SNESNGMRESFormCombinedSolution_Private(snes,ivec,l,XM,FM,fMnorm,X,XA,FA);CHKERRQ(ierr);
262     /* r = F(x) */
263     if (fminnorm > fMnorm) fminnorm = fMnorm;  /* the minimum norm is now of F^M */
264 
265     /* differences for selection and restart */
266     if (ngmres->restart_type == SNES_NGMRES_RESTART_DIFFERENCE || ngmres->select_type == SNES_NGMRES_SELECT_DIFFERENCE) {
267       ierr = SNESNGMRESNorms_Private(snes,l,X,F,XM,FM,XA,FA,D,&dnorm,&dminnorm,&xMnorm,NULL,&yMnorm,&xAnorm,&fAnorm,&yAnorm);CHKERRQ(ierr);
268     } else {
269       ierr = SNESNGMRESNorms_Private(snes,l,X,F,XM,FM,XA,FA,D,NULL,NULL,&xMnorm,NULL,&yMnorm,&xAnorm,&fAnorm,&yAnorm);CHKERRQ(ierr);
270     }
271     SNESCheckFunctionNorm(snes,fnorm);
272 
273     /* combination (additive) or selection (multiplicative) of the N-GMRES solution */
274     ierr          = SNESNGMRESSelect_Private(snes,k_restart,XM,FM,xMnorm,fMnorm,yMnorm,XA,FA,xAnorm,fAnorm,yAnorm,dnorm,fminnorm,dminnorm,X,F,Y,&xnorm,&fnorm,&ynorm);CHKERRQ(ierr);
275     selectRestart = PETSC_FALSE;
276 
277     if (ngmres->restart_type == SNES_NGMRES_RESTART_DIFFERENCE) {
278       ierr = SNESNGMRESSelectRestart_Private(snes,l,fMnorm,fAnorm,dnorm,fminnorm,dminnorm,&selectRestart);CHKERRQ(ierr);
279 
280       /* if the restart conditions persist for more than restart_it iterations, restart. */
281       if (selectRestart) restart_count++;
282       else restart_count = 0;
283     } else if (ngmres->restart_type == SNES_NGMRES_RESTART_PERIODIC) {
284       if (k_restart > ngmres->restart_periodic) {
285         if (ngmres->monitor) ierr = PetscViewerASCIIPrintf(ngmres->monitor,"periodic restart after %D iterations\n",k_restart);CHKERRQ(ierr);
286         restart_count = ngmres->restart_it;
287       }
288     }
289 
290     ivec = k_restart % ngmres->msize; /* replace the last used part of the subspace */
291 
292     /* restart after restart conditions have persisted for a fixed number of iterations */
293     if (restart_count >= ngmres->restart_it) {
294       if (ngmres->monitor) {
295         ierr = PetscViewerASCIIPrintf(ngmres->monitor,"Restarted at iteration %d\n",k_restart);CHKERRQ(ierr);
296       }
297       restart_count = 0;
298       k_restart     = 1;
299       l             = 1;
300       ivec          = 0;
301       /* q_{00} = nu */
302       ierr = SNESNGMRESUpdateSubspace_Private(snes,0,0,FM,fMnorm,XM);CHKERRQ(ierr);
303     } else {
304       /* select the current size of the subspace */
305       if (l < ngmres->msize) l++;
306       k_restart++;
307       /* place the current entry in the list of previous entries */
308       if (ngmres->candidate) {
309         if (fminnorm > fMnorm) fminnorm = fMnorm;
310         ierr = SNESNGMRESUpdateSubspace_Private(snes,ivec,l,FM,fMnorm,XM);CHKERRQ(ierr);
311       } else {
312         if (fminnorm > fnorm) fminnorm = fnorm;
313         ierr = SNESNGMRESUpdateSubspace_Private(snes,ivec,l,F,fnorm,X);CHKERRQ(ierr);
314       }
315     }
316 
317     ierr       = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr);
318     snes->iter = k;
319     snes->norm = fnorm;
320     ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr);
321     ierr = SNESLogConvergenceHistory(snes,snes->norm,snes->iter);CHKERRQ(ierr);
322     ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr);
323     ierr = (*snes->ops->converged)(snes,snes->iter,0,0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr);
324     if (snes->reason) PetscFunctionReturn(0);
325   }
326   snes->reason = SNES_DIVERGED_MAX_IT;
327   PetscFunctionReturn(0);
328 }
329 
330 /*@
331  SNESNGMRESSetRestartFmRise - Increase the restart count if the step x_M increases the residual F_M
332 
333   Input Parameters:
334   +  snes - the SNES context.
335   -  flg  - boolean value deciding whether to use the option or not
336 
337   Options Database:
338   + -snes_ngmres_restart_fm_rise - Increase the restart count if the step x_M increases the residual F_M
339 
340   Level: intermediate
341 
342   Notes:
343   If the proposed step x_M increases the residual F_M, it might be trying to get out of a stagnation area.
344   To help the solver do that, reset the Krylov subspace whenever F_M increases.
345 
346   This option must be used with SNES_NGMRES_RESTART_DIFFERENCE
347 
348   The default is FALSE.
349   .seealso: SNES_NGMRES_RESTART_DIFFERENCE
350   @*/
351 PetscErrorCode SNESNGMRESSetRestartFmRise(SNES snes,PetscBool flg)
352 {
353     PetscErrorCode (*f)(SNES,PetscBool);
354     PetscErrorCode ierr;
355 
356     PetscFunctionBegin;
357     ierr = PetscObjectQueryFunction((PetscObject)snes,"SNESNGMRESSetRestartFmRise_C",&f);CHKERRQ(ierr);
358     if (f) {ierr = (f)(snes,flg);CHKERRQ(ierr);}
359     PetscFunctionReturn(0);
360 }
361 
362 PetscErrorCode SNESNGMRESSetRestartFmRise_NGMRES(SNES snes,PetscBool flg)
363 {
364   SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data;
365 
366   PetscFunctionBegin;
367   ngmres->restart_fm_rise = flg;
368   PetscFunctionReturn(0);
369 }
370 
371 PetscErrorCode SNESNGMRESGetRestartFmRise(SNES snes,PetscBool *flg)
372 {
373     PetscErrorCode (*f)(SNES,PetscBool*);
374     PetscErrorCode ierr;
375 
376     PetscFunctionBegin;
377     ierr = PetscObjectQueryFunction((PetscObject)snes,"SNESNGMRESGetRestartFmRise_C",&f);CHKERRQ(ierr);
378     if (f) {ierr = (f)(snes,flg);CHKERRQ(ierr);}
379     PetscFunctionReturn(0);
380 }
381 
382 PetscErrorCode SNESNGMRESGetRestartFmRise_NGMRES(SNES snes,PetscBool *flg)
383 {
384   SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data;
385 
386   PetscFunctionBegin;
387   *flg = ngmres->restart_fm_rise;
388   PetscFunctionReturn(0);
389 }
390 
391 
392 /*@
393     SNESNGMRESSetRestartType - Sets the restart type for SNESNGMRES.
394 
395     Logically Collective on SNES
396 
397     Input Parameters:
398 +   snes - the iterative context
399 -   rtype - restart type
400 
401     Options Database:
402 +   -snes_ngmres_restart_type<difference,periodic,none> - set the restart type
403 -   -snes_ngmres_restart[30] - sets the number of iterations before restart for periodic
404 
405     Level: intermediate
406 
407     SNESNGMRESRestartTypes:
408 +   SNES_NGMRES_RESTART_NONE - never restart
409 .   SNES_NGMRES_RESTART_DIFFERENCE - restart based upon difference criteria
410 -   SNES_NGMRES_RESTART_PERIODIC - restart after a fixed number of iterations
411 
412     Notes:
413     The default line search used is the L2 line search and it requires two additional function evaluations.
414 
415 @*/
416 PetscErrorCode SNESNGMRESSetRestartType(SNES snes,SNESNGMRESRestartType rtype)
417 {
418   PetscErrorCode ierr;
419 
420   PetscFunctionBegin;
421   PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
422   ierr = PetscTryMethod(snes,"SNESNGMRESSetRestartType_C",(SNES,SNESNGMRESRestartType),(snes,rtype));CHKERRQ(ierr);
423   PetscFunctionReturn(0);
424 }
425 
426 /*@
427     SNESNGMRESSetSelectType - Sets the selection type for SNESNGMRES.  This determines how the candidate solution and
428     combined solution are used to create the next iterate.
429 
430     Logically Collective on SNES
431 
432     Input Parameters:
433 +   snes - the iterative context
434 -   stype - selection type
435 
436     Options Database:
437 .   -snes_ngmres_select_type<difference,none,linesearch>
438 
439     Level: intermediate
440 
441     SNESNGMRESSelectTypes:
442 +   SNES_NGMRES_SELECT_NONE - choose the combined solution all the time
443 .   SNES_NGMRES_SELECT_DIFFERENCE - choose based upon the selection criteria
444 -   SNES_NGMRES_SELECT_LINESEARCH - choose based upon line search combination
445 
446     Notes:
447     The default line search used is the L2 line search and it requires two additional function evaluations.
448 
449 @*/
450 PetscErrorCode SNESNGMRESSetSelectType(SNES snes,SNESNGMRESSelectType stype)
451 {
452   PetscErrorCode ierr;
453 
454   PetscFunctionBegin;
455   PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
456   ierr = PetscTryMethod(snes,"SNESNGMRESSetSelectType_C",(SNES,SNESNGMRESSelectType),(snes,stype));CHKERRQ(ierr);
457   PetscFunctionReturn(0);
458 }
459 
460 PetscErrorCode SNESNGMRESSetSelectType_NGMRES(SNES snes,SNESNGMRESSelectType stype)
461 {
462   SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data;
463 
464   PetscFunctionBegin;
465   ngmres->select_type = stype;
466   PetscFunctionReturn(0);
467 }
468 
469 PetscErrorCode SNESNGMRESSetRestartType_NGMRES(SNES snes,SNESNGMRESRestartType rtype)
470 {
471   SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data;
472 
473   PetscFunctionBegin;
474   ngmres->restart_type = rtype;
475   PetscFunctionReturn(0);
476 }
477 
478 /*MC
479   SNESNGMRES - The Nonlinear Generalized Minimum Residual method.
480 
481    Level: beginner
482 
483    Options Database:
484 +  -snes_ngmres_select_type<difference,none,linesearch> - choose the select between candidate and combined solution
485 .  -snes_ngmres_restart_type<difference,none,periodic> - choose the restart conditions
486 .  -snes_ngmres_candidate        - Use NGMRES variant which combines candidate solutions instead of actual solutions
487 .  -snes_ngmres_m                - Number of stored previous solutions and residuals
488 .  -snes_ngmres_restart_it       - Number of iterations the restart conditions hold before restart
489 .  -snes_ngmres_gammaA           - Residual tolerance for solution select between the candidate and combination
490 .  -snes_ngmres_gammaC           - Residual tolerance for restart
491 .  -snes_ngmres_epsilonB         - Difference tolerance between subsequent solutions triggering restart
492 .  -snes_ngmres_deltaB           - Difference tolerance between residuals triggering restart
493 .  -snes_ngmres_restart_fm_rise  - Restart on residual rise from x_M step
494 .  -snes_ngmres_monitor          - Prints relevant information about the ngmres iteration
495 .  -snes_linesearch_type <basic,l2,cp> - Line search type used for the default smoother
496 -  -additive_snes_linesearch_type - linesearch type used to select between the candidate and combined solution with additive select type
497 
498    Notes:
499 
500    The N-GMRES method combines m previous solutions into a minimum-residual solution by solving a small linearized
501    optimization problem at each iteration.
502 
503    Very similar to the SNESANDERSON algorithm.
504 
505    References:
506 +  1. - C. W. Oosterlee and T. Washio, "Krylov Subspace Acceleration of Nonlinear Multigrid with Application to Recirculating Flows",
507    SIAM Journal on Scientific Computing, 21(5), 2000.
508 -  2. - Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu, "Composing Scalable Nonlinear Algebraic Solvers",
509    SIAM Review, 57(4), 2015
510 
511 
512 .seealso: SNESCreate(), SNES, SNESSetType(), SNESType (for list of available types)
513 M*/
514 
515 PETSC_EXTERN PetscErrorCode SNESCreate_NGMRES(SNES snes)
516 {
517   SNES_NGMRES    *ngmres;
518   PetscErrorCode ierr;
519   SNESLineSearch linesearch;
520 
521   PetscFunctionBegin;
522   snes->ops->destroy        = SNESDestroy_NGMRES;
523   snes->ops->setup          = SNESSetUp_NGMRES;
524   snes->ops->setfromoptions = SNESSetFromOptions_NGMRES;
525   snes->ops->view           = SNESView_NGMRES;
526   snes->ops->solve          = SNESSolve_NGMRES;
527   snes->ops->reset          = SNESReset_NGMRES;
528 
529   snes->usesnpc  = PETSC_TRUE;
530   snes->usesksp  = PETSC_FALSE;
531   snes->npcside  = PC_RIGHT;
532 
533   snes->alwayscomputesfinalresidual = PETSC_TRUE;
534 
535   ierr          = PetscNewLog(snes,&ngmres);CHKERRQ(ierr);
536   snes->data    = (void*) ngmres;
537   ngmres->msize = 30;
538 
539   if (!snes->tolerancesset) {
540     snes->max_funcs = 30000;
541     snes->max_its   = 10000;
542   }
543 
544   ngmres->candidate = PETSC_FALSE;
545 
546   ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr);
547   if (!((PetscObject)linesearch)->type_name) {
548     ierr = SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC);CHKERRQ(ierr);
549   }
550 
551   ngmres->additive_linesearch = NULL;
552   ngmres->approxfunc          = PETSC_FALSE;
553   ngmres->restart_it          = 2;
554   ngmres->restart_periodic    = 30;
555   ngmres->gammaA              = 2.0;
556   ngmres->gammaC              = 2.0;
557   ngmres->deltaB              = 0.9;
558   ngmres->epsilonB            = 0.1;
559   ngmres->restart_fm_rise     = PETSC_FALSE;
560 
561   ngmres->restart_type = SNES_NGMRES_RESTART_DIFFERENCE;
562   ngmres->select_type  = SNES_NGMRES_SELECT_DIFFERENCE;
563 
564   ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetSelectType_C",SNESNGMRESSetSelectType_NGMRES);CHKERRQ(ierr);
565   ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetRestartType_C",SNESNGMRESSetRestartType_NGMRES);CHKERRQ(ierr);
566   ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetRestartFmRise_C",SNESNGMRESSetRestartFmRise_NGMRES);CHKERRQ(ierr);
567   ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESGetRestartFmRise_C",SNESNGMRESGetRestartFmRise_NGMRES);CHKERRQ(ierr);
568   PetscFunctionReturn(0);
569 }
570