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