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