xref: /petsc/src/snes/impls/ngmres/snesngmres.c (revision 72edecb9bd3aeb4a9296b816611f1c42fbe134ed)
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   PCSide              npcside;
184   PetscErrorCode      ierr;
185 
186   PetscFunctionBegin;
187   /* variable initialization */
188   snes->reason = SNES_CONVERGED_ITERATING;
189   X            = snes->vec_sol;
190   F            = snes->vec_func;
191   B            = snes->vec_rhs;
192   XA           = snes->vec_sol_update;
193   FA           = snes->work[0];
194   D            = snes->work[1];
195 
196   /* work for the line search */
197   Y  = snes->work[2];
198   XM = snes->work[3];
199   FM = snes->work[4];
200 
201   ierr       = PetscObjectAMSTakeAccess((PetscObject)snes);CHKERRQ(ierr);
202   snes->iter = 0;
203   snes->norm = 0.;
204   ierr       = PetscObjectAMSGrantAccess((PetscObject)snes);CHKERRQ(ierr);
205 
206   /* initialization */
207 
208   if (snes->pc && snes->pcside == PC_LEFT) {
209     ierr = SNESApplyPC(snes,X,NULL,NULL,F);CHKERRQ(ierr);
210     ierr = SNESGetConvergedReason(snes->pc,&reason);CHKERRQ(ierr);
211     if (reason < 0  && reason != SNES_DIVERGED_MAX_IT) {
212       snes->reason = SNES_DIVERGED_INNER;
213       PetscFunctionReturn(0);
214     }
215     ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr);
216   } else {
217     if (!snes->vec_func_init_set) {
218       ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr);
219       if (snes->domainerror) {
220         snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN;
221         PetscFunctionReturn(0);
222       }
223     } else snes->vec_func_init_set = PETSC_FALSE;
224     if (!snes->norm_init_set) {
225       /* convergence test */
226       ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr);
227       if (PetscIsInfOrNanReal(fnorm)) {
228         snes->reason = SNES_DIVERGED_FNORM_NAN;
229         PetscFunctionReturn(0);
230       }
231     } else {
232       fnorm               = snes->norm_init;
233       snes->norm_init_set = PETSC_FALSE;
234     }
235   }
236   fminnorm = fnorm;
237 
238   /* q_{00} = nu  */
239   SNESNGMRESUpdateSubspace_Private(snes,0,0,F,fnorm,X);
240 
241   ierr       = PetscObjectAMSTakeAccess((PetscObject)snes);CHKERRQ(ierr);
242   snes->norm = fnorm;
243   ierr       = PetscObjectAMSGrantAccess((PetscObject)snes);CHKERRQ(ierr);
244   ierr       = SNESLogConvergenceHistory(snes,fnorm,0);CHKERRQ(ierr);
245   ierr       = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr);
246   ierr       = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr);
247   if (snes->reason) PetscFunctionReturn(0);
248 
249   k_restart = 1;
250   l         = 1;
251   for (k=1; k < snes->max_its+1; k++) {
252     /* select which vector of the stored subspace will be updated */
253     ivec = k_restart % ngmres->msize; /* replace the last used part of the subspace */
254 
255     /* Computation of x^M */
256     if (snes->pc && snes->pcside == PC_RIGHT) {
257       ierr = VecCopy(X,XM);CHKERRQ(ierr);
258       ierr = SNESSetInitialFunction(snes->pc,F);CHKERRQ(ierr);
259       ierr = SNESSetInitialFunctionNorm(snes->pc,fnorm);CHKERRQ(ierr);
260 
261       ierr = PetscLogEventBegin(SNES_NPCSolve,snes->pc,XM,B,0);CHKERRQ(ierr);
262       ierr = SNESSolve(snes->pc,B,XM);CHKERRQ(ierr);
263       ierr = PetscLogEventEnd(SNES_NPCSolve,snes->pc,XM,B,0);CHKERRQ(ierr);
264 
265       ierr = SNESGetConvergedReason(snes->pc,&reason);CHKERRQ(ierr);
266       if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
267         snes->reason = SNES_DIVERGED_INNER;
268         PetscFunctionReturn(0);
269       }
270       ierr = SNESGetPCFunction(snes,FM,&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 
509   ierr          = PetscNewLog(snes,SNES_NGMRES,&ngmres);CHKERRQ(ierr);
510   snes->data    = (void*) ngmres;
511   ngmres->msize = 30;
512 
513   if (!snes->tolerancesset) {
514     snes->max_funcs = 30000;
515     snes->max_its   = 10000;
516   }
517 
518   ngmres->candidate = PETSC_FALSE;
519 
520   ngmres->additive_linesearch = NULL;
521   ngmres->approxfunc       = PETSC_FALSE;
522   ngmres->restart_it       = 2;
523   ngmres->restart_periodic = 30;
524   ngmres->gammaA           = 2.0;
525   ngmres->gammaC           = 2.0;
526   ngmres->deltaB           = 0.9;
527   ngmres->epsilonB         = 0.1;
528 
529   ngmres->restart_type = SNES_NGMRES_RESTART_DIFFERENCE;
530   ngmres->select_type  = SNES_NGMRES_SELECT_DIFFERENCE;
531 
532   ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetSelectType_C",SNESNGMRESSetSelectType_NGMRES);CHKERRQ(ierr);
533   ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetRestartType_C",SNESNGMRESSetRestartType_NGMRES);CHKERRQ(ierr);
534   PetscFunctionReturn(0);
535 }
536 
537