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