xref: /petsc/src/snes/impls/ncg/snesncg.c (revision 5a856986583887c326abe5dfd149e8184a29cd80)
1 #include <../src/snes/impls/ncg/snesncgimpl.h> /*I "petscsnes.h" I*/
2 const char *const SNESNCGTypes[] = {"FR","PRP","HS","DY","CD","SNESNCGType","SNES_NCG_",0};
3 
4 PetscErrorCode SNESReset_NCG(SNES snes)
5 {
6   PetscFunctionBegin;
7   PetscFunctionReturn(0);
8 }
9 
10 #define SNESLINESEARCHNCGLINEAR "ncglinear"
11 
12 /*
13   SNESDestroy_NCG - Destroys the private SNES_NCG context that was created with SNESCreate_NCG().
14 
15   Input Parameter:
16 . snes - the SNES context
17 
18   Application Interface Routine: SNESDestroy()
19 */
20 PetscErrorCode SNESDestroy_NCG(SNES snes)
21 {
22   PetscErrorCode ierr;
23 
24   PetscFunctionBegin;
25   ierr = PetscFree(snes->data);CHKERRQ(ierr);
26   PetscFunctionReturn(0);
27 }
28 
29 /*
30    SNESSetUp_NCG - Sets up the internal data structures for the later use
31    of the SNESNCG nonlinear solver.
32 
33    Input Parameters:
34 +  snes - the SNES context
35 -  x - the solution vector
36 
37    Application Interface Routine: SNESSetUp()
38  */
39 
40 PETSC_EXTERN PetscErrorCode SNESLineSearchCreate_NCGLinear(SNESLineSearch);
41 
42 PetscErrorCode SNESSetUp_NCG(SNES snes)
43 {
44   PetscErrorCode ierr;
45 
46   PetscFunctionBegin;
47   ierr = SNESSetWorkVecs(snes,2);CHKERRQ(ierr);
48   if (snes->npcside== PC_RIGHT) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "SNESNCG only supports left preconditioning");
49   if (snes->functype == SNES_FUNCTION_DEFAULT) snes->functype = SNES_FUNCTION_UNPRECONDITIONED;
50   PetscFunctionReturn(0);
51 }
52 /*
53   SNESSetFromOptions_NCG - Sets various parameters for the SNESNCG method.
54 
55   Input Parameter:
56 . snes - the SNES context
57 
58   Application Interface Routine: SNESSetFromOptions()
59 */
60 static PetscErrorCode SNESSetFromOptions_NCG(PetscOptionItems *PetscOptionsObject,SNES snes)
61 {
62   SNES_NCG       *ncg = (SNES_NCG*)snes->data;
63   PetscErrorCode ierr;
64   PetscBool      debug = PETSC_FALSE;
65   SNESNCGType    ncgtype=ncg->type;
66   SNESLineSearch linesearch;
67 
68   PetscFunctionBegin;
69   ierr = PetscOptionsHead(PetscOptionsObject,"SNES NCG options");CHKERRQ(ierr);
70   ierr = PetscOptionsBool("-snes_ncg_monitor","Monitor NCG iterations","SNES",ncg->monitor ? PETSC_TRUE : PETSC_FALSE, &debug, NULL);CHKERRQ(ierr);
71   if (debug) {
72     ncg->monitor = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes));CHKERRQ(ierr);
73   }
74   ierr = PetscOptionsEnum("-snes_ncg_type","NCG Beta type used","SNESNCGSetType",SNESNCGTypes,(PetscEnum)ncg->type,(PetscEnum*)&ncgtype,NULL);CHKERRQ(ierr);
75   ierr = SNESNCGSetType(snes, ncgtype);CHKERRQ(ierr);
76   ierr = PetscOptionsTail();CHKERRQ(ierr);
77   if (!snes->linesearch) {
78     ierr = SNESGetLineSearch(snes, &linesearch);CHKERRQ(ierr);
79     if (!snes->npc) {
80       ierr = SNESLineSearchSetType(linesearch, SNESLINESEARCHCP);CHKERRQ(ierr);
81     } else {
82       ierr = SNESLineSearchSetType(linesearch, SNESLINESEARCHL2);CHKERRQ(ierr);
83     }
84   }
85   ierr = SNESLineSearchRegister(SNESLINESEARCHNCGLINEAR, SNESLineSearchCreate_NCGLinear);CHKERRQ(ierr);
86   PetscFunctionReturn(0);
87 }
88 
89 /*
90   SNESView_NCG - Prints info from the SNESNCG data structure.
91 
92   Input Parameters:
93 + SNES - the SNES context
94 - viewer - visualization context
95 
96   Application Interface Routine: SNESView()
97 */
98 static PetscErrorCode SNESView_NCG(SNES snes, PetscViewer viewer)
99 {
100   SNES_NCG      *ncg = (SNES_NCG *) snes->data;
101   PetscBool      iascii;
102   PetscErrorCode ierr;
103 
104   PetscFunctionBegin;
105   ierr = PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr);
106   if (iascii) {
107     ierr = PetscViewerASCIIPrintf(viewer, "  type: %s\n", SNESNCGTypes[ncg->type]);CHKERRQ(ierr);
108   }
109   PetscFunctionReturn(0);
110 }
111 
112 PetscErrorCode SNESLineSearchApply_NCGLinear(SNESLineSearch linesearch)
113 {
114   PetscScalar    alpha, ptAp;
115   Vec            X, Y, F, W;
116   SNES           snes;
117   PetscErrorCode ierr;
118   PetscReal      *fnorm, *xnorm, *ynorm;
119 
120   PetscFunctionBegin;
121   ierr  = SNESLineSearchGetSNES(linesearch, &snes);CHKERRQ(ierr);
122   X     = linesearch->vec_sol;
123   W     = linesearch->vec_sol_new;
124   F     = linesearch->vec_func;
125   Y     = linesearch->vec_update;
126   fnorm = &linesearch->fnorm;
127   xnorm = &linesearch->xnorm;
128   ynorm = &linesearch->ynorm;
129 
130   if (!snes->jacobian) {
131     ierr = SNESSetUpMatrices(snes);CHKERRQ(ierr);
132   }
133 
134   /*
135 
136    The exact step size for unpreconditioned linear CG is just:
137    alpha = (r, r) / (p, Ap) = (f, f) / (y, Jy)
138    */
139   ierr  = SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre);CHKERRQ(ierr);
140   ierr  = VecDot(F, F, &alpha);CHKERRQ(ierr);
141   ierr  = MatMult(snes->jacobian, Y, W);CHKERRQ(ierr);
142   ierr  = VecDot(Y, W, &ptAp);CHKERRQ(ierr);
143   alpha = alpha / ptAp;
144   ierr  = VecAXPY(X,-alpha,Y);CHKERRQ(ierr);
145   ierr  = SNESComputeFunction(snes,X,F);CHKERRQ(ierr);
146 
147   ierr = VecNorm(F,NORM_2,fnorm);CHKERRQ(ierr);
148   ierr = VecNorm(X,NORM_2,xnorm);CHKERRQ(ierr);
149   ierr = VecNorm(Y,NORM_2,ynorm);CHKERRQ(ierr);
150   PetscFunctionReturn(0);
151 }
152 
153 PETSC_EXTERN PetscErrorCode SNESLineSearchCreate_NCGLinear(SNESLineSearch linesearch)
154 {
155   PetscFunctionBegin;
156   linesearch->ops->apply          = SNESLineSearchApply_NCGLinear;
157   linesearch->ops->destroy        = NULL;
158   linesearch->ops->setfromoptions = NULL;
159   linesearch->ops->reset          = NULL;
160   linesearch->ops->view           = NULL;
161   linesearch->ops->setup          = NULL;
162   PetscFunctionReturn(0);
163 }
164 
165 /*
166 
167  Assuming F = SNESComputeFunction(X) compute Y^tJ^tF using a simple secant approximation of the jacobian.
168 
169  */
170 PetscErrorCode SNESNCGComputeYtJtF_Private(SNES snes, Vec X, Vec F, Vec Y, Vec W, Vec G, PetscScalar * ytJtf)
171 {
172   PetscErrorCode ierr;
173   PetscScalar    ftf, ftg, fty, h;
174 
175   PetscFunctionBegin;
176   ierr   = VecDot(F, F, &ftf);CHKERRQ(ierr);
177   ierr   = VecDot(F, Y, &fty);CHKERRQ(ierr);
178   h      = 1e-5*fty / fty;
179   ierr   = VecCopy(X, W);CHKERRQ(ierr);
180   ierr   = VecAXPY(W, -h, Y);CHKERRQ(ierr);          /* this is arbitrary */
181   ierr   = SNESComputeFunction(snes, W, G);CHKERRQ(ierr);
182   ierr   = VecDot(G, F, &ftg);CHKERRQ(ierr);
183   *ytJtf = (ftg - ftf) / h;
184   PetscFunctionReturn(0);
185 }
186 
187 /*@
188     SNESNCGSetType - Sets the conjugate update type for SNESNCG.
189 
190     Logically Collective on SNES
191 
192     Input Parameters:
193 +   snes - the iterative context
194 -   btype - update type
195 
196     Options Database:
197 .   -snes_ncg_type<prp,fr,hs,dy,cd>
198 
199     Level: intermediate
200 
201     SNESNCGSelectTypes:
202 +   SNES_NCG_FR - Fletcher-Reeves update
203 .   SNES_NCG_PRP - Polak-Ribiere-Polyak update
204 .   SNES_NCG_HS - Hestenes-Steifel update
205 .   SNES_NCG_DY - Dai-Yuan update
206 -   SNES_NCG_CD - Conjugate Descent update
207 
208    Notes:
209    PRP is the default, and the only one that tolerates generalized search directions.
210 
211 @*/
212 PetscErrorCode SNESNCGSetType(SNES snes, SNESNCGType btype)
213 {
214   PetscErrorCode ierr;
215 
216   PetscFunctionBegin;
217   PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
218   ierr = PetscTryMethod(snes,"SNESNCGSetType_C",(SNES,SNESNCGType),(snes,btype));CHKERRQ(ierr);
219   PetscFunctionReturn(0);
220 }
221 
222 PetscErrorCode SNESNCGSetType_NCG(SNES snes, SNESNCGType btype)
223 {
224   SNES_NCG *ncg = (SNES_NCG*)snes->data;
225 
226   PetscFunctionBegin;
227   ncg->type = btype;
228   PetscFunctionReturn(0);
229 }
230 
231 /*
232   SNESSolve_NCG - Solves a nonlinear system with the Nonlinear Conjugate Gradient method.
233 
234   Input Parameters:
235 . snes - the SNES context
236 
237   Output Parameter:
238 . outits - number of iterations until termination
239 
240   Application Interface Routine: SNESSolve()
241 */
242 PetscErrorCode SNESSolve_NCG(SNES snes)
243 {
244   SNES_NCG             *ncg = (SNES_NCG*)snes->data;
245   Vec                  X,dX,lX,F,dXold;
246   PetscReal            fnorm, ynorm, xnorm, beta = 0.0;
247   PetscScalar          dXdotdX, dXolddotdXold, dXdotdXold, lXdotdX, lXdotdXold;
248   PetscInt             maxits, i;
249   PetscErrorCode       ierr;
250   SNESLineSearchReason lsresult = SNES_LINESEARCH_SUCCEEDED;
251   SNESLineSearch       linesearch;
252   SNESConvergedReason  reason;
253 
254   PetscFunctionBegin;
255   if (snes->xl || snes->xu || snes->ops->computevariablebounds) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name);
256 
257   ierr = PetscCitationsRegister(SNESCitation,&SNEScite);CHKERRQ(ierr);
258   snes->reason = SNES_CONVERGED_ITERATING;
259 
260   maxits = snes->max_its;            /* maximum number of iterations */
261   X      = snes->vec_sol;            /* X^n */
262   dXold  = snes->work[0];            /* The previous iterate of X */
263   dX     = snes->work[1];            /* the preconditioned direction */
264   lX     = snes->vec_sol_update;     /* the conjugate direction */
265   F      = snes->vec_func;           /* residual vector */
266 
267   ierr = SNESGetLineSearch(snes, &linesearch);CHKERRQ(ierr);
268 
269   ierr       = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr);
270   snes->iter = 0;
271   snes->norm = 0.;
272   ierr       = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr);
273 
274   /* compute the initial function and preconditioned update dX */
275 
276   if (snes->npc && snes->functype == SNES_FUNCTION_PRECONDITIONED) {
277     ierr = SNESApplyNPC(snes,X,NULL,dX);CHKERRQ(ierr);
278     ierr = SNESGetConvergedReason(snes->npc,&reason);CHKERRQ(ierr);
279     if (reason < 0  && reason != SNES_DIVERGED_MAX_IT) {
280       snes->reason = SNES_DIVERGED_INNER;
281       PetscFunctionReturn(0);
282     }
283     ierr = VecCopy(dX,F);CHKERRQ(ierr);
284     ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr);
285   } else {
286     if (!snes->vec_func_init_set) {
287       ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr);
288     } else snes->vec_func_init_set = PETSC_FALSE;
289 
290     /* convergence test */
291     ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr);
292     SNESCheckFunctionNorm(snes,fnorm);
293     ierr = VecCopy(F,dX);CHKERRQ(ierr);
294   }
295   if (snes->npc) {
296     if (snes->functype == SNES_FUNCTION_UNPRECONDITIONED) {
297       ierr = SNESApplyNPC(snes,X,F,dX);CHKERRQ(ierr);
298       ierr = SNESGetConvergedReason(snes->npc,&reason);CHKERRQ(ierr);
299       if (reason < 0  && reason != SNES_DIVERGED_MAX_IT) {
300         snes->reason = SNES_DIVERGED_INNER;
301         PetscFunctionReturn(0);
302       }
303     }
304   }
305   ierr = VecCopy(dX,lX);CHKERRQ(ierr);
306   ierr = VecDot(dX, dX, &dXdotdX);CHKERRQ(ierr);
307 
308 
309   ierr       = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr);
310   snes->norm = fnorm;
311   ierr       = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr);
312   ierr       = SNESLogConvergenceHistory(snes,fnorm,0);CHKERRQ(ierr);
313   ierr       = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr);
314 
315   /* test convergence */
316   ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr);
317   if (snes->reason) PetscFunctionReturn(0);
318 
319   /* Call general purpose update function */
320   if (snes->ops->update) {
321     ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr);
322   }
323 
324   /* first update -- just use the (preconditioned) residual direction for the initial conjugate direction */
325 
326   for (i = 1; i < maxits + 1; i++) {
327     /* some update types require the old update direction or conjugate direction */
328     if (ncg->type != SNES_NCG_FR) {
329       ierr = VecCopy(dX, dXold);CHKERRQ(ierr);
330     }
331     ierr = SNESLineSearchApply(linesearch,X,F,&fnorm,lX);CHKERRQ(ierr);
332     ierr = SNESLineSearchGetReason(linesearch, &lsresult);CHKERRQ(ierr);
333     ierr = SNESLineSearchGetNorms(linesearch, &xnorm, &fnorm, &ynorm);CHKERRQ(ierr);
334     if (lsresult) {
335       if (++snes->numFailures >= snes->maxFailures) {
336         snes->reason = SNES_DIVERGED_LINE_SEARCH;
337         PetscFunctionReturn(0);
338       }
339     }
340     if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) {
341       snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
342       PetscFunctionReturn(0);
343     }
344     /* Monitor convergence */
345     ierr       = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr);
346     snes->iter = i;
347     snes->norm = fnorm;
348     snes->xnorm = xnorm;
349     snes->ynorm = ynorm;
350     ierr       = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr);
351     ierr       = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr);
352     ierr       = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr);
353 
354     /* Test for convergence */
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     /* Call general purpose update function */
359     if (snes->ops->update) {
360       ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr);
361     }
362     if (snes->npc) {
363       if (snes->functype == SNES_FUNCTION_PRECONDITIONED) {
364         ierr = SNESApplyNPC(snes,X,NULL,dX);CHKERRQ(ierr);
365         ierr = SNESGetConvergedReason(snes->npc,&reason);CHKERRQ(ierr);
366         if (reason < 0  && reason != SNES_DIVERGED_MAX_IT) {
367           snes->reason = SNES_DIVERGED_INNER;
368           PetscFunctionReturn(0);
369         }
370         ierr = VecCopy(dX,F);CHKERRQ(ierr);
371       } else {
372         ierr = SNESApplyNPC(snes,X,F,dX);CHKERRQ(ierr);
373         ierr = SNESGetConvergedReason(snes->npc,&reason);CHKERRQ(ierr);
374         if (reason < 0  && reason != SNES_DIVERGED_MAX_IT) {
375           snes->reason = SNES_DIVERGED_INNER;
376           PetscFunctionReturn(0);
377         }
378       }
379     } else {
380       ierr = VecCopy(F,dX);CHKERRQ(ierr);
381     }
382 
383     /* compute the conjugate direction lX = dX + beta*lX with beta = ((dX, dX) / (dX_old, dX_old) (Fletcher-Reeves update)*/
384     switch (ncg->type) {
385     case SNES_NCG_FR: /* Fletcher-Reeves */
386       dXolddotdXold= dXdotdX;
387       ierr         = VecDot(dX, dX, &dXdotdX);CHKERRQ(ierr);
388       beta         = PetscRealPart(dXdotdX / dXolddotdXold);
389       break;
390     case SNES_NCG_PRP: /* Polak-Ribiere-Poylak */
391       dXolddotdXold= dXdotdX;
392       ierr         = VecDotBegin(dX, dX, &dXdotdXold);CHKERRQ(ierr);
393       ierr         = VecDotBegin(dXold, dX, &dXdotdXold);CHKERRQ(ierr);
394       ierr         = VecDotEnd(dX, dX, &dXdotdX);CHKERRQ(ierr);
395       ierr         = VecDotEnd(dXold, dX, &dXdotdXold);CHKERRQ(ierr);
396       beta         = PetscRealPart(((dXdotdX - dXdotdXold) / dXolddotdXold));
397       if (beta < 0.0) beta = 0.0; /* restart */
398       break;
399     case SNES_NCG_HS: /* Hestenes-Stiefel */
400       ierr = VecDotBegin(dX, dX, &dXdotdX);CHKERRQ(ierr);
401       ierr = VecDotBegin(dX, dXold, &dXdotdXold);CHKERRQ(ierr);
402       ierr = VecDotBegin(lX, dX, &lXdotdX);CHKERRQ(ierr);
403       ierr = VecDotBegin(lX, dXold, &lXdotdXold);CHKERRQ(ierr);
404       ierr = VecDotEnd(dX, dX, &dXdotdX);CHKERRQ(ierr);
405       ierr = VecDotEnd(dX, dXold, &dXdotdXold);CHKERRQ(ierr);
406       ierr = VecDotEnd(lX, dX, &lXdotdX);CHKERRQ(ierr);
407       ierr = VecDotEnd(lX, dXold, &lXdotdXold);CHKERRQ(ierr);
408       beta = PetscRealPart((dXdotdX - dXdotdXold) / (lXdotdX - lXdotdXold));
409       break;
410     case SNES_NCG_DY: /* Dai-Yuan */
411       ierr = VecDotBegin(dX, dX, &dXdotdX);CHKERRQ(ierr);
412       ierr = VecDotBegin(lX, dX, &lXdotdX);CHKERRQ(ierr);
413       ierr = VecDotBegin(lX, dXold, &lXdotdXold);CHKERRQ(ierr);
414       ierr = VecDotEnd(dX, dX, &dXdotdX);CHKERRQ(ierr);
415       ierr = VecDotEnd(lX, dX, &lXdotdX);CHKERRQ(ierr);
416       ierr = VecDotEnd(lX, dXold, &lXdotdXold);CHKERRQ(ierr);
417       beta = PetscRealPart(dXdotdX / (lXdotdXold - lXdotdX));CHKERRQ(ierr);
418       break;
419     case SNES_NCG_CD: /* Conjugate Descent */
420       ierr = VecDotBegin(dX, dX, &dXdotdX);CHKERRQ(ierr);
421       ierr = VecDotBegin(lX, dXold, &lXdotdXold);CHKERRQ(ierr);
422       ierr = VecDotEnd(dX, dX, &dXdotdX);CHKERRQ(ierr);
423       ierr = VecDotEnd(lX, dXold, &lXdotdXold);CHKERRQ(ierr);
424       beta = PetscRealPart(dXdotdX / lXdotdXold);CHKERRQ(ierr);
425       break;
426     }
427     if (ncg->monitor) {
428       ierr = PetscViewerASCIIPrintf(ncg->monitor, "beta = %e\n", (double)beta);CHKERRQ(ierr);
429     }
430     ierr = VecAYPX(lX, beta, dX);CHKERRQ(ierr);
431   }
432   ierr = PetscInfo1(snes, "Maximum number of iterations has been reached: %D\n", maxits);CHKERRQ(ierr);
433   if (!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT;
434   PetscFunctionReturn(0);
435 }
436 
437 
438 
439 /*MC
440   SNESNCG - Nonlinear Conjugate-Gradient method for the solution of general nonlinear systems.
441 
442   Level: beginner
443 
444   Options Database:
445 +   -snes_ncg_type <fr, prp, dy, hs, cd> - Choice of conjugate-gradient update parameter, default is prp.
446 .   -snes_linesearch_type <cp,l2,basic> - Line search type.
447 -   -snes_ncg_monitor - Print relevant information about the ncg iteration.
448 
449    Notes:
450     This solves the nonlinear system of equations F(x) = 0 using the nonlinear generalization of the conjugate
451           gradient method.  This may be used with a nonlinear preconditioner used to pick the new search directions, but otherwise
452           chooses the initial search direction as F(x) for the initial guess x.
453 
454           Only supports left non-linear preconditioning.
455 
456    References:
457 .  1. -  Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu,"Composing Scalable Nonlinear Algebraic Solvers",
458    SIAM Review, 57(4), 2015
459 
460 
461 .seealso:  SNESCreate(), SNES, SNESSetType(), SNESNEWTONLS, SNESNEWTONTR, SNESNGMRES, SNESQN
462 M*/
463 PETSC_EXTERN PetscErrorCode SNESCreate_NCG(SNES snes)
464 {
465   PetscErrorCode ierr;
466   SNES_NCG       * neP;
467 
468   PetscFunctionBegin;
469   snes->ops->destroy        = SNESDestroy_NCG;
470   snes->ops->setup          = SNESSetUp_NCG;
471   snes->ops->setfromoptions = SNESSetFromOptions_NCG;
472   snes->ops->view           = SNESView_NCG;
473   snes->ops->solve          = SNESSolve_NCG;
474   snes->ops->reset          = SNESReset_NCG;
475 
476   snes->usesksp = PETSC_FALSE;
477   snes->usesnpc = PETSC_TRUE;
478   snes->npcside = PC_LEFT;
479 
480   snes->alwayscomputesfinalresidual = PETSC_TRUE;
481 
482   if (!snes->tolerancesset) {
483     snes->max_funcs = 30000;
484     snes->max_its   = 10000;
485     snes->stol      = 1e-20;
486   }
487 
488   ierr         = PetscNewLog(snes,&neP);CHKERRQ(ierr);
489   snes->data   = (void*) neP;
490   neP->monitor = NULL;
491   neP->type    = SNES_NCG_PRP;
492   ierr         = PetscObjectComposeFunction((PetscObject)snes,"SNESNCGSetType_C", SNESNCGSetType_NCG);CHKERRQ(ierr);
493   PetscFunctionReturn(0);
494 }
495