Actual source code: ls.c
1: #define PETSCSNES_DLL
3: #include ../src/snes/impls/ls/ls.h
5: /*
6: Checks if J^T F = 0 which implies we've found a local minimum of the norm of the function,
7: || F(u) ||_2 but not a zero, F(u) = 0. In the case when one cannot compute J^T F we use the fact that
8: 0 = (J^T F)^T W = F^T J W iff W not in the null space of J. Thanks for Jorge More
9: for this trick. One assumes that the probability that W is in the null space of J is very, very small.
10: */
13: PetscErrorCode SNESLSCheckLocalMin_Private(SNES snes,Mat A,Vec F,Vec W,PetscReal fnorm,PetscTruth *ismin)
14: {
15: PetscReal a1;
17: PetscTruth hastranspose;
20: *ismin = PETSC_FALSE;
21: MatHasOperation(A,MATOP_MULT_TRANSPOSE,&hastranspose);
22: if (hastranspose) {
23: /* Compute || J^T F|| */
24: MatMultTranspose(A,F,W);
25: VecNorm(W,NORM_2,&a1);
26: PetscInfo1(snes,"|| J^T F|| %G near zero implies found a local minimum\n",a1/fnorm);
27: if (a1/fnorm < 1.e-4) *ismin = PETSC_TRUE;
28: } else {
29: Vec work;
30: PetscScalar result;
31: PetscReal wnorm;
33: VecSetRandom(W,PETSC_NULL);
34: VecNorm(W,NORM_2,&wnorm);
35: VecDuplicate(W,&work);
36: MatMult(A,W,work);
37: VecDot(F,work,&result);
38: VecDestroy(work);
39: a1 = PetscAbsScalar(result)/(fnorm*wnorm);
40: PetscInfo1(snes,"(F^T J random)/(|| F ||*||random|| %G near zero implies found a local minimum\n",a1);
41: if (a1 < 1.e-4) *ismin = PETSC_TRUE;
42: }
43: return(0);
44: }
46: /*
47: Checks if J^T(F - J*X) = 0
48: */
51: PetscErrorCode SNESLSCheckResidual_Private(SNES snes,Mat A,Vec F,Vec X,Vec W1,Vec W2)
52: {
53: PetscReal a1,a2;
55: PetscTruth hastranspose;
58: MatHasOperation(A,MATOP_MULT_TRANSPOSE,&hastranspose);
59: if (hastranspose) {
60: MatMult(A,X,W1);
61: VecAXPY(W1,-1.0,F);
63: /* Compute || J^T W|| */
64: MatMultTranspose(A,W1,W2);
65: VecNorm(W1,NORM_2,&a1);
66: VecNorm(W2,NORM_2,&a2);
67: if (a1 != 0.0) {
68: PetscInfo1(snes,"||J^T(F-Ax)||/||F-AX|| %G near zero implies inconsistent rhs\n",a2/a1);
69: }
70: }
71: return(0);
72: }
74: /* --------------------------------------------------------------------
76: This file implements a truncated Newton method with a line search,
77: for solving a system of nonlinear equations, using the KSP, Vec,
78: and Mat interfaces for linear solvers, vectors, and matrices,
79: respectively.
81: The following basic routines are required for each nonlinear solver:
82: SNESCreate_XXX() - Creates a nonlinear solver context
83: SNESSetFromOptions_XXX() - Sets runtime options
84: SNESSolve_XXX() - Solves the nonlinear system
85: SNESDestroy_XXX() - Destroys the nonlinear solver context
86: The suffix "_XXX" denotes a particular implementation, in this case
87: we use _LS (e.g., SNESCreate_LS, SNESSolve_LS) for solving
88: systems of nonlinear equations with a line search (LS) method.
89: These routines are actually called via the common user interface
90: routines SNESCreate(), SNESSetFromOptions(), SNESSolve(), and
91: SNESDestroy(), so the application code interface remains identical
92: for all nonlinear solvers.
94: Another key routine is:
95: SNESSetUp_XXX() - Prepares for the use of a nonlinear solver
96: by setting data structures and options. The interface routine SNESSetUp()
97: is not usually called directly by the user, but instead is called by
98: SNESSolve() if necessary.
100: Additional basic routines are:
101: SNESView_XXX() - Prints details of runtime options that
102: have actually been used.
103: These are called by application codes via the interface routines
104: SNESView().
106: The various types of solvers (preconditioners, Krylov subspace methods,
107: nonlinear solvers, timesteppers) are all organized similarly, so the
108: above description applies to these categories also.
110: -------------------------------------------------------------------- */
111: /*
112: SNESSolve_LS - Solves a nonlinear system with a truncated Newton
113: method with a line search.
115: Input Parameters:
116: . snes - the SNES context
118: Output Parameter:
119: . outits - number of iterations until termination
121: Application Interface Routine: SNESSolve()
123: Notes:
124: This implements essentially a truncated Newton method with a
125: line search. By default a cubic backtracking line search
126: is employed, as described in the text "Numerical Methods for
127: Unconstrained Optimization and Nonlinear Equations" by Dennis
128: and Schnabel.
129: */
132: PetscErrorCode SNESSolve_LS(SNES snes)
133: {
134: SNES_LS *neP = (SNES_LS*)snes->data;
135: PetscErrorCode ierr;
136: PetscInt maxits,i,lits;
137: PetscTruth lssucceed;
138: MatStructure flg = DIFFERENT_NONZERO_PATTERN;
139: PetscReal fnorm,gnorm,xnorm=0,ynorm;
140: Vec Y,X,F,G,W;
141: KSPConvergedReason kspreason;
144: snes->numFailures = 0;
145: snes->numLinearSolveFailures = 0;
146: snes->reason = SNES_CONVERGED_ITERATING;
148: maxits = snes->max_its; /* maximum number of iterations */
149: X = snes->vec_sol; /* solution vector */
150: F = snes->vec_func; /* residual vector */
151: Y = snes->work[0]; /* work vectors */
152: G = snes->work[1];
153: W = snes->work[2];
155: PetscObjectTakeAccess(snes);
156: snes->iter = 0;
157: snes->norm = 0.0;
158: PetscObjectGrantAccess(snes);
159: SNESComputeFunction(snes,X,F);
160: if (snes->domainerror) {
161: snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN;
162: return(0);
163: }
164: VecNormBegin(F,NORM_2,&fnorm); /* fnorm <- ||F|| */
165: VecNormBegin(X,NORM_2,&xnorm); /* xnorm <- ||x|| */
166: VecNormEnd(F,NORM_2,&fnorm);
167: VecNormEnd(X,NORM_2,&xnorm);
168: if PetscIsInfOrNanReal(fnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
169: PetscObjectTakeAccess(snes);
170: snes->norm = fnorm;
171: PetscObjectGrantAccess(snes);
172: SNESLogConvHistory(snes,fnorm,0);
173: SNESMonitor(snes,0,fnorm);
175: /* set parameter for default relative tolerance convergence test */
176: snes->ttol = fnorm*snes->rtol;
177: /* test convergence */
178: (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);
179: if (snes->reason) return(0);
181: for (i=0; i<maxits; i++) {
183: /* Call general purpose update function */
184: if (snes->ops->update) {
185: (*snes->ops->update)(snes, snes->iter);
186: }
188: /* Solve J Y = F, where J is Jacobian matrix */
189: SNESComputeJacobian(snes,X,&snes->jacobian,&snes->jacobian_pre,&flg);
190: KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre,flg);
191: SNES_KSPSolve(snes,snes->ksp,F,Y);
192: KSPGetConvergedReason(snes->ksp,&kspreason);
193: if (kspreason < 0) {
194: if (++snes->numLinearSolveFailures >= snes->maxLinearSolveFailures) {
195: PetscInfo2(snes,"iter=%D, number linear solve failures %D greater than current SNES allowed, stopping solve\n",snes->iter,snes->numLinearSolveFailures);
196: snes->reason = SNES_DIVERGED_LINEAR_SOLVE;
197: break;
198: }
199: }
200: KSPGetIterationNumber(snes->ksp,&lits);
201: snes->linear_its += lits;
202: PetscInfo2(snes,"iter=%D, linear solve iterations=%D\n",snes->iter,lits);
204: if (neP->precheckstep) {
205: PetscTruth changed_y = PETSC_FALSE;
206: (*neP->precheckstep)(snes,X,Y,neP->precheck,&changed_y);
207: }
209: if (PetscLogPrintInfo){
210: SNESLSCheckResidual_Private(snes,snes->jacobian,F,Y,G,W);
211: }
213: /* Compute a (scaled) negative update in the line search routine:
214: Y <- X - lambda*Y
215: and evaluate G = function(Y) (depends on the line search).
216: */
217: VecCopy(Y,snes->vec_sol_update);
218: ynorm = 1; gnorm = fnorm;
219: (*neP->LineSearch)(snes,neP->lsP,X,F,G,Y,W,fnorm,xnorm,&ynorm,&gnorm,&lssucceed);
220: PetscInfo4(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lssucceed=%d\n",fnorm,gnorm,ynorm,(int)lssucceed);
221: if (snes->reason == SNES_DIVERGED_FUNCTION_COUNT) break;
222: if (snes->domainerror) {
223: snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN;
224: return(0);
225: }
226: if (!lssucceed) {
227: if (++snes->numFailures >= snes->maxFailures) {
228: PetscTruth ismin;
229: snes->reason = SNES_DIVERGED_LS_FAILURE;
230: SNESLSCheckLocalMin_Private(snes,snes->jacobian,G,W,gnorm,&ismin);
231: if (ismin) snes->reason = SNES_DIVERGED_LOCAL_MIN;
232: break;
233: }
234: }
235: /* Update function and solution vectors */
236: fnorm = gnorm;
237: VecCopy(G,F);
238: VecCopy(W,X);
239: /* Monitor convergence */
240: PetscObjectTakeAccess(snes);
241: snes->iter = i+1;
242: snes->norm = fnorm;
243: PetscObjectGrantAccess(snes);
244: SNESLogConvHistory(snes,snes->norm,lits);
245: SNESMonitor(snes,snes->iter,snes->norm);
246: /* Test for convergence, xnorm = || X || */
247: if (snes->ops->converged != SNESSkipConverged) { VecNorm(X,NORM_2,&xnorm); }
248: (*snes->ops->converged)(snes,snes->iter,xnorm,ynorm,fnorm,&snes->reason,snes->cnvP);
249: if (snes->reason) break;
250: }
251: if (i == maxits) {
252: PetscInfo1(snes,"Maximum number of iterations has been reached: %D\n",maxits);
253: if(!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT;
254: }
255: return(0);
256: }
257: /* -------------------------------------------------------------------------- */
258: /*
259: SNESSetUp_LS - Sets up the internal data structures for the later use
260: of the SNESLS nonlinear solver.
262: Input Parameter:
263: . snes - the SNES context
264: . x - the solution vector
266: Application Interface Routine: SNESSetUp()
268: Notes:
269: For basic use of the SNES solvers, the user need not explicitly call
270: SNESSetUp(), since these actions will automatically occur during
271: the call to SNESSolve().
272: */
275: PetscErrorCode SNESSetUp_LS(SNES snes)
276: {
280: if (!snes->vec_sol_update) {
281: VecDuplicate(snes->vec_sol,&snes->vec_sol_update);
282: PetscLogObjectParent(snes,snes->vec_sol_update);
283: }
284: if (!snes->work) {
285: snes->nwork = 3;
286: VecDuplicateVecs(snes->vec_sol,snes->nwork,&snes->work);
287: PetscLogObjectParents(snes,snes->nwork,snes->work);
288: }
289: return(0);
290: }
291: /* -------------------------------------------------------------------------- */
292: /*
293: SNESDestroy_LS - Destroys the private SNES_LS context that was created
294: with SNESCreate_LS().
296: Input Parameter:
297: . snes - the SNES context
299: Application Interface Routine: SNESDestroy()
300: */
303: PetscErrorCode SNESDestroy_LS(SNES snes)
304: {
308: if (snes->vec_sol_update) {
309: VecDestroy(snes->vec_sol_update);
310: snes->vec_sol_update = PETSC_NULL;
311: }
312: if (snes->nwork) {
313: VecDestroyVecs(snes->work,snes->nwork);
314: snes->nwork = 0;
315: snes->work = PETSC_NULL;
316: }
317: PetscFree(snes->data);
319: /* clear composed functions */
320: PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSet_C","",PETSC_NULL);
321: PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSetPostCheck_C","",PETSC_NULL);
322: PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSetPreCheck_C","",PETSC_NULL);
324: return(0);
325: }
326: /* -------------------------------------------------------------------------- */
330: /*@C
331: SNESLineSearchNo - This routine is not a line search at all;
332: it simply uses the full Newton step. Thus, this routine is intended
333: to serve as a template and is not recommended for general use.
335: Collective on SNES and Vec
337: Input Parameters:
338: + snes - nonlinear context
339: . lsctx - optional context for line search (not used here)
340: . x - current iterate
341: . f - residual evaluated at x
342: . y - search direction
343: . fnorm - 2-norm of f
344: - xnorm - norm of x if known, otherwise 0
346: Output Parameters:
347: + g - residual evaluated at new iterate y
348: . w - new iterate
349: . gnorm - 2-norm of g
350: . ynorm - 2-norm of search length
351: - flag - PETSC_TRUE on success, PETSC_FALSE on failure
353: Options Database Key:
354: . -snes_ls basic - Activates SNESLineSearchNo()
356: Level: advanced
358: .keywords: SNES, nonlinear, line search, cubic
360: .seealso: SNESLineSearchCubic(), SNESLineSearchQuadratic(),
361: SNESLineSearchSet(), SNESLineSearchNoNorms()
362: @*/
363: PetscErrorCode SNESLineSearchNo(SNES snes,void *lsctx,Vec x,Vec f,Vec g,Vec y,Vec w,PetscReal fnorm,PetscReal xnorm,PetscReal *ynorm,PetscReal *gnorm,PetscTruth *flag)
364: {
366: SNES_LS *neP = (SNES_LS*)snes->data;
367: PetscTruth changed_w = PETSC_FALSE,changed_y = PETSC_FALSE;
370: *flag = PETSC_TRUE;
371: PetscLogEventBegin(SNES_LineSearch,snes,x,f,g);
372: VecNorm(y,NORM_2,ynorm); /* ynorm = || y || */
373: VecWAXPY(w,-1.0,y,x); /* w <- x - y */
374: if (neP->postcheckstep) {
375: (*neP->postcheckstep)(snes,x,y,w,neP->postcheck,&changed_y,&changed_w);
376: }
377: if (changed_y) {
378: VecWAXPY(w,-1.0,y,x); /* w <- x - y */
379: }
380: SNESComputeFunction(snes,w,g);
381: if (!snes->domainerror) {
382: VecNorm(g,NORM_2,gnorm); /* gnorm = || g || */
383: if PetscIsInfOrNanReal(*gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
384: }
385: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
386: return(0);
387: }
388: /* -------------------------------------------------------------------------- */
392: /*@C
393: SNESLineSearchNoNorms - This routine is not a line search at
394: all; it simply uses the full Newton step. This version does not
395: even compute the norm of the function or search direction; this
396: is intended only when you know the full step is fine and are
397: not checking for convergence of the nonlinear iteration (for
398: example, you are running always for a fixed number of Newton steps).
400: Collective on SNES and Vec
402: Input Parameters:
403: + snes - nonlinear context
404: . lsctx - optional context for line search (not used here)
405: . x - current iterate
406: . f - residual evaluated at x
407: . y - search direction
408: . w - work vector
409: . fnorm - 2-norm of f
410: - xnorm - norm of x if known, otherwise 0
412: Output Parameters:
413: + g - residual evaluated at new iterate y
414: . w - new iterate
415: . gnorm - not changed
416: . ynorm - not changed
417: - flag - set to PETSC_TRUE indicating a successful line search
419: Options Database Key:
420: . -snes_ls basicnonorms - Activates SNESLineSearchNoNorms()
422: Notes:
423: SNESLineSearchNoNorms() must be used in conjunction with
424: either the options
425: $ -snes_no_convergence_test -snes_max_it <its>
426: or alternatively a user-defined custom test set via
427: SNESSetConvergenceTest(); or a -snes_max_it of 1,
428: otherwise, the SNES solver will generate an error.
430: During the final iteration this will not evaluate the function at
431: the solution point. This is to save a function evaluation while
432: using pseudo-timestepping.
434: The residual norms printed by monitoring routines such as
435: SNESMonitorDefault() (as activated via -snes_monitor) will not be
436: correct, since they are not computed.
438: Level: advanced
440: .keywords: SNES, nonlinear, line search, cubic
442: .seealso: SNESLineSearchCubic(), SNESLineSearchQuadratic(),
443: SNESLineSearchSet(), SNESLineSearchNo()
444: @*/
445: PetscErrorCode SNESLineSearchNoNorms(SNES snes,void *lsctx,Vec x,Vec f,Vec g,Vec y,Vec w,PetscReal fnorm,PetscReal xnorm,PetscReal *ynorm,PetscReal *gnorm,PetscTruth *flag)
446: {
448: SNES_LS *neP = (SNES_LS*)snes->data;
449: PetscTruth changed_w = PETSC_FALSE,changed_y = PETSC_FALSE;
452: *flag = PETSC_TRUE;
453: PetscLogEventBegin(SNES_LineSearch,snes,x,f,g);
454: VecWAXPY(w,-1.0,y,x); /* w <- x - y */
455: if (neP->postcheckstep) {
456: (*neP->postcheckstep)(snes,x,y,w,neP->postcheck,&changed_y,&changed_w);
457: }
458: if (changed_y) {
459: VecWAXPY(w,-1.0,y,x); /* w <- x - y */
460: }
461:
462: /* don't evaluate function the last time through */
463: if (snes->iter < snes->max_its-1) {
464: SNESComputeFunction(snes,w,g);
465: }
466: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
467: return(0);
468: }
469: /* -------------------------------------------------------------------------- */
472: /*@C
473: SNESLineSearchCubic - Performs a cubic line search (default line search method).
475: Collective on SNES
477: Input Parameters:
478: + snes - nonlinear context
479: . lsctx - optional context for line search (not used here)
480: . x - current iterate
481: . f - residual evaluated at x
482: . y - search direction
483: . w - work vector
484: . fnorm - 2-norm of f
485: - xnorm - norm of x if known, otherwise 0
487: Output Parameters:
488: + g - residual evaluated at new iterate y
489: . w - new iterate
490: . gnorm - 2-norm of g
491: . ynorm - 2-norm of search length
492: - flag - PETSC_TRUE if line search succeeds; PETSC_FALSE on failure.
494: Options Database Key:
495: + -snes_ls cubic - Activates SNESLineSearchCubic()
496: . -snes_ls_alpha <alpha> - Sets alpha
497: . -snes_ls_maxstep <maxstep> - Sets the maximum stepsize the line search will use (if the 2-norm(y) > maxstep then scale y to be y = (maxstep/2-norm(y)) *y)
498: - -snes_ls_minlambda <minlambda> - Sets the minimum lambda the line search will use minlambda/ max_i ( y[i]/x[i] )
500:
501: Notes:
502: This line search is taken from "Numerical Methods for Unconstrained
503: Optimization and Nonlinear Equations" by Dennis and Schnabel, page 325.
505: Level: advanced
507: .keywords: SNES, nonlinear, line search, cubic
509: .seealso: SNESLineSearchQuadratic(), SNESLineSearchNo(), SNESLineSearchSet(), SNESLineSearchNoNorms()
510: @*/
511: PetscErrorCode SNESLineSearchCubic(SNES snes,void *lsctx,Vec x,Vec f,Vec g,Vec y,Vec w,PetscReal fnorm,PetscReal xnorm,PetscReal *ynorm,PetscReal *gnorm,PetscTruth *flag)
512: {
513: /*
514: Note that for line search purposes we work with with the related
515: minimization problem:
516: min z(x): R^n -> R,
517: where z(x) = .5 * fnorm*fnorm, and fnorm = || f ||_2.
518: */
519:
520: PetscReal initslope,lambdaprev,gnormprev,a,b,d,t1,t2,rellength;
521: PetscReal minlambda,lambda,lambdatemp;
522: #if defined(PETSC_USE_COMPLEX)
523: PetscScalar cinitslope;
524: #endif
526: PetscInt count;
527: SNES_LS *neP = (SNES_LS*)snes->data;
528: PetscTruth changed_w = PETSC_FALSE,changed_y = PETSC_FALSE;
531: PetscLogEventBegin(SNES_LineSearch,snes,x,f,g);
532: *flag = PETSC_TRUE;
534: VecNorm(y,NORM_2,ynorm);
535: if (*ynorm == 0.0) {
536: PetscInfo(snes,"Search direction and size is 0\n");
537: *gnorm = fnorm;
538: VecCopy(x,w);
539: VecCopy(f,g);
540: *flag = PETSC_FALSE;
541: goto theend1;
542: }
543: if (*ynorm > neP->maxstep) { /* Step too big, so scale back */
544: PetscInfo2(snes,"Scaling step by %G old ynorm %G\n",neP->maxstep/(*ynorm),*ynorm);
545: VecScale(y,neP->maxstep/(*ynorm));
546: *ynorm = neP->maxstep;
547: }
548: VecMaxPointwiseDivide(y,x,&rellength);
549: minlambda = neP->minlambda/rellength;
550: MatMult(snes->jacobian,y,w);
551: #if defined(PETSC_USE_COMPLEX)
552: VecDot(f,w,&cinitslope);
553: initslope = PetscRealPart(cinitslope);
554: #else
555: VecDot(f,w,&initslope);
556: #endif
557: if (initslope > 0.0) initslope = -initslope;
558: if (initslope == 0.0) initslope = -1.0;
560: VecWAXPY(w,-1.0,y,x);
561: if (snes->nfuncs >= snes->max_funcs) {
562: PetscInfo(snes,"Exceeded maximum function evaluations, while checking full step length!\n");
563: *flag = PETSC_FALSE;
564: snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
565: goto theend1;
566: }
567: SNESComputeFunction(snes,w,g);
568: if (snes->domainerror) {
569: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
570: return(0);
571: }
572: VecNorm(g,NORM_2,gnorm);
573: if PetscIsInfOrNanReal(*gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
574: PetscInfo2(snes,"Initial fnorm %G gnorm %G\n",fnorm,*gnorm);
575: if (.5*(*gnorm)*(*gnorm) <= .5*fnorm*fnorm + neP->alpha*initslope) { /* Sufficient reduction */
576: PetscInfo2(snes,"Using full step: fnorm %G gnorm %G\n",fnorm,*gnorm);
577: goto theend1;
578: }
580: /* Fit points with quadratic */
581: lambda = 1.0;
582: lambdatemp = -initslope/((*gnorm)*(*gnorm) - fnorm*fnorm - 2.0*initslope);
583: lambdaprev = lambda;
584: gnormprev = *gnorm;
585: if (lambdatemp > .5*lambda) lambdatemp = .5*lambda;
586: if (lambdatemp <= .1*lambda) lambda = .1*lambda;
587: else lambda = lambdatemp;
589: VecWAXPY(w,-lambda,y,x);
590: if (snes->nfuncs >= snes->max_funcs) {
591: PetscInfo1(snes,"Exceeded maximum function evaluations, while attempting quadratic backtracking! %D \n",snes->nfuncs);
592: *flag = PETSC_FALSE;
593: snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
594: goto theend1;
595: }
596: SNESComputeFunction(snes,w,g);
597: if (snes->domainerror) {
598: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
599: return(0);
600: }
601: VecNorm(g,NORM_2,gnorm);
602: if PetscIsInfOrNanReal(*gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
603: PetscInfo1(snes,"gnorm after quadratic fit %G\n",*gnorm);
604: if (.5*(*gnorm)*(*gnorm) < .5*fnorm*fnorm + lambda*neP->alpha*initslope) { /* sufficient reduction */
605: PetscInfo1(snes,"Quadratically determined step, lambda=%18.16e\n",lambda);
606: goto theend1;
607: }
609: /* Fit points with cubic */
610: count = 1;
611: while (PETSC_TRUE) {
612: if (lambda <= minlambda) {
613: PetscInfo1(snes,"Unable to find good step length! After %D tries \n",count);
614: PetscInfo6(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, minlambda=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",fnorm,*gnorm,*ynorm,minlambda,lambda,initslope);
615: *flag = PETSC_FALSE;
616: break;
617: }
618: t1 = .5*((*gnorm)*(*gnorm) - fnorm*fnorm) - lambda*initslope;
619: t2 = .5*(gnormprev*gnormprev - fnorm*fnorm) - lambdaprev*initslope;
620: a = (t1/(lambda*lambda) - t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev);
621: b = (-lambdaprev*t1/(lambda*lambda) + lambda*t2/(lambdaprev*lambdaprev))/(lambda-lambdaprev);
622: d = b*b - 3*a*initslope;
623: if (d < 0.0) d = 0.0;
624: if (a == 0.0) {
625: lambdatemp = -initslope/(2.0*b);
626: } else {
627: lambdatemp = (-b + sqrt(d))/(3.0*a);
628: }
629: lambdaprev = lambda;
630: gnormprev = *gnorm;
631: if (lambdatemp > .5*lambda) lambdatemp = .5*lambda;
632: if (lambdatemp <= .1*lambda) lambda = .1*lambda;
633: else lambda = lambdatemp;
634: VecWAXPY(w,-lambda,y,x);
635: if (snes->nfuncs >= snes->max_funcs) {
636: PetscInfo1(snes,"Exceeded maximum function evaluations, while looking for good step length! %D \n",count);
637: PetscInfo5(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",fnorm,*gnorm,*ynorm,lambda,initslope);
638: *flag = PETSC_FALSE;
639: snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
640: break;
641: }
642: SNESComputeFunction(snes,w,g);
643: if (snes->domainerror) {
644: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
645: return(0);
646: }
647: VecNorm(g,NORM_2,gnorm);
648: if PetscIsInfOrNanReal(*gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
649: if (.5*(*gnorm)*(*gnorm) < .5*fnorm*fnorm + lambda*neP->alpha*initslope) { /* is reduction enough? */
650: PetscInfo2(snes,"Cubically determined step, current gnorm %G lambda=%18.16e\n",*gnorm,lambda);
651: break;
652: } else {
653: PetscInfo2(snes,"Cubic step no good, shrinking lambda, current gnorem %G lambda=%18.16e\n",*gnorm,lambda);
654: }
655: count++;
656: }
657: theend1:
658: /* Optional user-defined check for line search step validity */
659: if (neP->postcheckstep && *flag) {
660: (*neP->postcheckstep)(snes,x,y,w,neP->postcheck,&changed_y,&changed_w);
661: if (changed_y) {
662: VecWAXPY(w,-1.0,y,x);
663: }
664: if (changed_y || changed_w) { /* recompute the function if the step has changed */
665: SNESComputeFunction(snes,w,g);
666: if (snes->domainerror) {
667: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
668: return(0);
669: }
670: VecNormBegin(g,NORM_2,gnorm);
671: if PetscIsInfOrNanReal(*gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
672: VecNormBegin(y,NORM_2,ynorm);
673: VecNormEnd(g,NORM_2,gnorm);
674: VecNormEnd(y,NORM_2,ynorm);
675: }
676: }
677: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
678: return(0);
679: }
680: /* -------------------------------------------------------------------------- */
683: /*@C
684: SNESLineSearchQuadratic - Performs a quadratic line search.
686: Collective on SNES and Vec
688: Input Parameters:
689: + snes - the SNES context
690: . lsctx - optional context for line search (not used here)
691: . x - current iterate
692: . f - residual evaluated at x
693: . y - search direction
694: . w - work vector
695: . fnorm - 2-norm of f
696: - xnorm - norm of x if known, otherwise 0
698: Output Parameters:
699: + g - residual evaluated at new iterate w
700: . w - new iterate (x + lambda*y)
701: . gnorm - 2-norm of g
702: . ynorm - 2-norm of search length
703: - flag - PETSC_TRUE if line search succeeds; PETSC_FALSE on failure.
705: Options Database Keys:
706: + -snes_ls quadratic - Activates SNESLineSearchQuadratic()
707: . -snes_ls_alpha <alpha> - Sets alpha
708: . -snes_ls_maxstep <maxstep> - Sets the maximum stepsize the line search will use (if the 2-norm(y) > maxstep then scale y to be y = (maxstep/2-norm(y)) *y)
709: - -snes_ls_minlambda <minlambda> - Sets the minimum lambda the line search will use minlambda/ max_i ( y[i]/x[i] )
711: Notes:
712: Use SNESLineSearchSet() to set this routine within the SNESLS method.
714: Level: advanced
716: .keywords: SNES, nonlinear, quadratic, line search
718: .seealso: SNESLineSearchCubic(), SNESLineSearchNo(), SNESLineSearchSet(), SNESLineSearchNoNorms()
719: @*/
720: PetscErrorCode SNESLineSearchQuadratic(SNES snes,void *lsctx,Vec x,Vec f,Vec g,Vec y,Vec w,PetscReal fnorm,PetscReal xnorm,PetscReal *ynorm,PetscReal *gnorm,PetscTruth *flag)
721: {
722: /*
723: Note that for line search purposes we work with with the related
724: minimization problem:
725: min z(x): R^n -> R,
726: where z(x) = .5 * fnorm*fnorm,and fnorm = || f ||_2.
727: */
728: PetscReal initslope,minlambda,lambda,lambdatemp,rellength;
729: #if defined(PETSC_USE_COMPLEX)
730: PetscScalar cinitslope;
731: #endif
733: PetscInt count;
734: SNES_LS *neP = (SNES_LS*)snes->data;
735: PetscTruth changed_w = PETSC_FALSE,changed_y = PETSC_FALSE;
738: PetscLogEventBegin(SNES_LineSearch,snes,x,f,g);
739: *flag = PETSC_TRUE;
741: VecNorm(y,NORM_2,ynorm);
742: if (*ynorm == 0.0) {
743: PetscInfo(snes,"Search direction and size is 0\n");
744: *gnorm = fnorm;
745: VecCopy(x,w);
746: VecCopy(f,g);
747: *flag = PETSC_FALSE;
748: goto theend2;
749: }
750: if (*ynorm > neP->maxstep) { /* Step too big, so scale back */
751: VecScale(y,neP->maxstep/(*ynorm));
752: *ynorm = neP->maxstep;
753: }
754: VecMaxPointwiseDivide(y,x,&rellength);
755: minlambda = neP->minlambda/rellength;
756: MatMult(snes->jacobian,y,w);
757: #if defined(PETSC_USE_COMPLEX)
758: VecDot(f,w,&cinitslope);
759: initslope = PetscRealPart(cinitslope);
760: #else
761: VecDot(f,w,&initslope);
762: #endif
763: if (initslope > 0.0) initslope = -initslope;
764: if (initslope == 0.0) initslope = -1.0;
765: PetscInfo1(snes,"Initslope %G \n",initslope);
767: VecWAXPY(w,-1.0,y,x);
768: if (snes->nfuncs >= snes->max_funcs) {
769: PetscInfo(snes,"Exceeded maximum function evaluations, while checking full step length!\n");
770: *flag = PETSC_FALSE;
771: snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
772: goto theend2;
773: }
774: SNESComputeFunction(snes,w,g);
775: if (snes->domainerror) {
776: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
777: return(0);
778: }
779: VecNorm(g,NORM_2,gnorm);
780: if PetscIsInfOrNanReal(*gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
781: if (.5*(*gnorm)*(*gnorm) <= .5*fnorm*fnorm + neP->alpha*initslope) { /* Sufficient reduction */
782: PetscInfo(snes,"Using full step\n");
783: goto theend2;
784: }
786: /* Fit points with quadratic */
787: lambda = 1.0;
788: count = 1;
789: while (PETSC_TRUE) {
790: if (lambda <= minlambda) { /* bad luck; use full step */
791: PetscInfo1(snes,"Unable to find good step length! %D \n",count);
792: PetscInfo5(snes,"fnorm=%G, gnorm=%G, ynorm=%G, lambda=%G, initial slope=%G\n",fnorm,*gnorm,*ynorm,lambda,initslope);
793: VecCopy(x,w);
794: *flag = PETSC_FALSE;
795: break;
796: }
797: lambdatemp = -initslope/((*gnorm)*(*gnorm) - fnorm*fnorm - 2.0*initslope);
798: if (lambdatemp > .5*lambda) lambdatemp = .5*lambda;
799: if (lambdatemp <= .1*lambda) lambda = .1*lambda;
800: else lambda = lambdatemp;
801:
802: VecWAXPY(w,-lambda,y,x);
803: if (snes->nfuncs >= snes->max_funcs) {
804: PetscInfo1(snes,"Exceeded maximum function evaluations, while looking for good step length! %D \n",count);
805: PetscInfo5(snes,"fnorm=%18.16e, gnorm=%18.16e, ynorm=%18.16e, lambda=%18.16e, initial slope=%18.16e\n",fnorm,*gnorm,*ynorm,lambda,initslope);
806: *flag = PETSC_FALSE;
807: snes->reason = SNES_DIVERGED_FUNCTION_COUNT;
808: break;
809: }
810: SNESComputeFunction(snes,w,g);
811: if (snes->domainerror) {
812: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
813: return(0);
814: }
815: VecNorm(g,NORM_2,gnorm);
816: if PetscIsInfOrNanReal(*gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
817: if (.5*(*gnorm)*(*gnorm) < .5*fnorm*fnorm + lambda*neP->alpha*initslope) { /* sufficient reduction */
818: PetscInfo1(snes,"Quadratically determined step, lambda=%G\n",lambda);
819: break;
820: }
821: count++;
822: }
823: theend2:
824: /* Optional user-defined check for line search step validity */
825: if (neP->postcheckstep) {
826: (*neP->postcheckstep)(snes,x,y,w,neP->postcheck,&changed_y,&changed_w);
827: if (changed_y) {
828: VecWAXPY(w,-1.0,y,x);
829: }
830: if (changed_y || changed_w) { /* recompute the function if the step has changed */
831: SNESComputeFunction(snes,w,g);
832: if (snes->domainerror) {
833: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
834: return(0);
835: }
836: VecNormBegin(g,NORM_2,gnorm);
837: VecNormBegin(y,NORM_2,ynorm);
838: VecNormEnd(g,NORM_2,gnorm);
839: VecNormEnd(y,NORM_2,ynorm);
840: if PetscIsInfOrNanReal(*gnorm) SETERRQ(PETSC_ERR_FP,"User provided compute function generated a Not-a-Number");
841: }
842: }
843: PetscLogEventEnd(SNES_LineSearch,snes,x,f,g);
844: return(0);
845: }
847: /* -------------------------------------------------------------------------- */
850: /*@C
851: SNESLineSearchSet - Sets the line search routine to be used
852: by the method SNESLS.
854: Input Parameters:
855: + snes - nonlinear context obtained from SNESCreate()
856: . lsctx - optional user-defined context for use by line search
857: - func - pointer to int function
859: Collective on SNES
861: Available Routines:
862: + SNESLineSearchCubic() - default line search
863: . SNESLineSearchQuadratic() - quadratic line search
864: . SNESLineSearchNo() - the full Newton step (actually not a line search)
865: - SNESLineSearchNoNorms() - the full Newton step (calculating no norms; faster in parallel)
867: Options Database Keys:
868: + -snes_ls [cubic,quadratic,basic,basicnonorms] - Selects line search
869: . -snes_ls_alpha <alpha> - Sets alpha
870: . -snes_ls_maxstep <maxstep> - Sets maximum step the line search will use (if the 2-norm(y) > maxstep then scale y to be y = (maxstep/2-norm(y)) *y)
871: - -snes_ls_minlambda <minlambda> - Sets the minimum lambda the line search will use minlambda / max_i ( y[i]/x[i] )
873: Calling sequence of func:
874: .vb
875: func (SNES snes,void *lsctx,Vec x,Vec f,Vec g,Vec y,Vec w,PetscReal fnorm,PetscReal xnorm,PetscReal *ynorm,PetscReal *gnorm,PetscTruth *flag)
876: .ve
878: Input parameters for func:
879: + snes - nonlinear context
880: . lsctx - optional user-defined context for line search
881: . x - current iterate
882: . f - residual evaluated at x
883: . y - search direction
884: - fnorm - 2-norm of f
886: Output parameters for func:
887: + g - residual evaluated at new iterate y
888: . w - new iterate
889: . gnorm - 2-norm of g
890: . ynorm - 2-norm of search length
891: - flag - set to PETSC_TRUE if the line search succeeds; PETSC_FALSE on failure.
893: Level: advanced
895: .keywords: SNES, nonlinear, set, line search, routine
897: .seealso: SNESLineSearchCubic(), SNESLineSearchQuadratic(), SNESLineSearchNo(), SNESLineSearchNoNorms(),
898: SNESLineSearchSetPostCheck(), SNESLineSearchSetParams(), SNESLineSearchGetParams(), SNESLineSearchSetPreCheck()
899: @*/
900: PetscErrorCode SNESLineSearchSet(SNES snes,PetscErrorCode (*func)(SNES,void*,Vec,Vec,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal*,PetscReal*,PetscTruth*),void *lsctx)
901: {
902: PetscErrorCode ierr,(*f)(SNES,PetscErrorCode (*)(SNES,void*,Vec,Vec,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal*,PetscReal*,PetscTruth*),void*);
905: PetscObjectQueryFunction((PetscObject)snes,"SNESLineSearchSet_C",(void (**)(void))&f);
906: if (f) {
907: (*f)(snes,func,lsctx);
908: }
909: return(0);
910: }
913: /* -------------------------------------------------------------------------- */
917: PetscErrorCode SNESLineSearchSet_LS(SNES snes,FCN2 func,void *lsctx)
918: {
920: ((SNES_LS *)(snes->data))->LineSearch = func;
921: ((SNES_LS *)(snes->data))->lsP = lsctx;
922: return(0);
923: }
925: /* -------------------------------------------------------------------------- */
928: /*@C
929: SNESLineSearchSetPostCheck - Sets a routine to check the validity of new iterate computed
930: by the line search routine in the Newton-based method SNESLS.
932: Input Parameters:
933: + snes - nonlinear context obtained from SNESCreate()
934: . func - pointer to function
935: - checkctx - optional user-defined context for use by step checking routine
937: Collective on SNES
939: Calling sequence of func:
940: .vb
941: int func (SNES snes, Vec x,Vec y,Vec w,void *checkctx, PetscTruth *changed_y,PetscTruth *changed_w)
942: .ve
943: where func returns an error code of 0 on success and a nonzero
944: on failure.
946: Input parameters for func:
947: + snes - nonlinear context
948: . checkctx - optional user-defined context for use by step checking routine
949: . x - previous iterate
950: . y - new search direction and length
951: - w - current candidate iterate
953: Output parameters for func:
954: + y - search direction (possibly changed)
955: . w - current iterate (possibly modified)
956: . changed_y - indicates search direction was changed by this routine
957: - changed_w - indicates current iterate was changed by this routine
959: Level: advanced
961: Notes: All line searches accept the new iterate computed by the line search checking routine.
963: Only one of changed_y and changed_w can be PETSC_TRUE
965: On input w = x + y
967: SNESLineSearchNo() and SNESLineSearchNoNorms() (1) compute a candidate iterate u_{i+1}, (2) pass control
968: to the checking routine, and then (3) compute the corresponding nonlinear
969: function f(u_{i+1}) with the (possibly altered) iterate u_{i+1}.
971: SNESLineSearchQuadratic() and SNESLineSearchCubic() (1) compute a candidate iterate u_{i+1} as well as a
972: candidate nonlinear function f(u_{i+1}), (2) pass control to the checking
973: routine, and then (3) force a re-evaluation of f(u_{i+1}) if any changes
974: were made to the candidate iterate in the checking routine (as indicated
975: by flag=PETSC_TRUE). The overhead of this extra function re-evaluation can be
976: very costly, so use this feature with caution!
978: .keywords: SNES, nonlinear, set, line search check, step check, routine
980: .seealso: SNESLineSearchSet(), SNESLineSearchSetPreCheck()
981: @*/
982: PetscErrorCode SNESLineSearchSetPostCheck(SNES snes,PetscErrorCode (*func)(SNES,Vec,Vec,Vec,void*,PetscTruth*,PetscTruth*),void *checkctx)
983: {
984: PetscErrorCode ierr,(*f)(SNES,PetscErrorCode (*)(SNES,Vec,Vec,Vec,void*,PetscTruth*,PetscTruth*),void*);
987: PetscObjectQueryFunction((PetscObject)snes,"SNESLineSearchSetPostCheck_C",(void (**)(void))&f);
988: if (f) {
989: (*f)(snes,func,checkctx);
990: }
991: return(0);
992: }
995: /*@C
996: SNESLineSearchSetPreCheck - Sets a routine to check the validity of a new direction given by the linear solve
997: before the line search is called.
999: Input Parameters:
1000: + snes - nonlinear context obtained from SNESCreate()
1001: . func - pointer to function
1002: - checkctx - optional user-defined context for use by step checking routine
1004: Collective on SNES
1006: Calling sequence of func:
1007: .vb
1008: int func (SNES snes, Vec x,Vec y,void *checkctx, PetscTruth *changed_y)
1009: .ve
1010: where func returns an error code of 0 on success and a nonzero
1011: on failure.
1013: Input parameters for func:
1014: + snes - nonlinear context
1015: . checkctx - optional user-defined context for use by step checking routine
1016: . x - previous iterate
1017: - y - new search direction and length
1019: Output parameters for func:
1020: + y - search direction (possibly changed)
1021: - changed_y - indicates search direction was changed by this routine
1023: Level: advanced
1025: Notes: All line searches accept the new iterate computed by the line search checking routine.
1027: .keywords: SNES, nonlinear, set, line search check, step check, routine
1029: .seealso: SNESLineSearchSet(), SNESLineSearchSetPostCheck(), SNESSetUpdate()
1030: @*/
1031: PetscErrorCode SNESLineSearchSetPreCheck(SNES snes,PetscErrorCode (*func)(SNES,Vec,Vec,void*,PetscTruth*),void *checkctx)
1032: {
1033: PetscErrorCode ierr,(*f)(SNES,PetscErrorCode (*)(SNES,Vec,Vec,void*,PetscTruth*),void*);
1036: PetscObjectQueryFunction((PetscObject)snes,"SNESLineSearchSetPreCheck_C",(void (**)(void))&f);
1037: if (f) {
1038: (*f)(snes,func,checkctx);
1039: }
1040: return(0);
1041: }
1043: /* -------------------------------------------------------------------------- */
1049: PetscErrorCode SNESLineSearchSetPostCheck_LS(SNES snes,FCN1 func,void *checkctx)
1050: {
1052: ((SNES_LS *)(snes->data))->postcheckstep = func;
1053: ((SNES_LS *)(snes->data))->postcheck = checkctx;
1054: return(0);
1055: }
1061: PetscErrorCode SNESLineSearchSetPreCheck_LS(SNES snes,FCN3 func,void *checkctx)
1062: {
1064: ((SNES_LS *)(snes->data))->precheckstep = func;
1065: ((SNES_LS *)(snes->data))->precheck = checkctx;
1066: return(0);
1067: }
1070: /*
1071: SNESView_LS - Prints info from the SNESLS data structure.
1073: Input Parameters:
1074: . SNES - the SNES context
1075: . viewer - visualization context
1077: Application Interface Routine: SNESView()
1078: */
1081: static PetscErrorCode SNESView_LS(SNES snes,PetscViewer viewer)
1082: {
1083: SNES_LS *ls = (SNES_LS *)snes->data;
1084: const char *cstr;
1086: PetscTruth iascii;
1089: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
1090: if (iascii) {
1091: if (ls->LineSearch == SNESLineSearchNo) cstr = "SNESLineSearchNo";
1092: else if (ls->LineSearch == SNESLineSearchQuadratic) cstr = "SNESLineSearchQuadratic";
1093: else if (ls->LineSearch == SNESLineSearchCubic) cstr = "SNESLineSearchCubic";
1094: else cstr = "unknown";
1095: PetscViewerASCIIPrintf(viewer," line search variant: %s\n",cstr);
1096: PetscViewerASCIIPrintf(viewer," alpha=%G, maxstep=%G, minlambda=%G\n",ls->alpha,ls->maxstep,ls->minlambda);
1097: } else {
1098: SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for SNES EQ LS",((PetscObject)viewer)->type_name);
1099: }
1100: return(0);
1101: }
1102: /* -------------------------------------------------------------------------- */
1103: /*
1104: SNESSetFromOptions_LS - Sets various parameters for the SNESLS method.
1106: Input Parameter:
1107: . snes - the SNES context
1109: Application Interface Routine: SNESSetFromOptions()
1110: */
1113: static PetscErrorCode SNESSetFromOptions_LS(SNES snes)
1114: {
1115: SNES_LS *ls = (SNES_LS *)snes->data;
1116: const char *lses[] = {"basic","basicnonorms","quadratic","cubic"};
1118: PetscInt indx;
1119: PetscTruth flg;
1122: PetscOptionsHead("SNES Line search options");
1123: PetscOptionsReal("-snes_ls_alpha","Function norm must decrease by","None",ls->alpha,&ls->alpha,0);
1124: PetscOptionsReal("-snes_ls_maxstep","Step must be less than","None",ls->maxstep,&ls->maxstep,0);
1125: PetscOptionsReal("-snes_ls_minlambda","Minimum lambda allowed","None",ls->minlambda,&ls->minlambda,0);
1127: PetscOptionsEList("-snes_ls","Line search used","SNESLineSearchSet",lses,4,"cubic",&indx,&flg);
1128: if (flg) {
1129: switch (indx) {
1130: case 0:
1131: SNESLineSearchSet(snes,SNESLineSearchNo,PETSC_NULL);
1132: break;
1133: case 1:
1134: SNESLineSearchSet(snes,SNESLineSearchNoNorms,PETSC_NULL);
1135: break;
1136: case 2:
1137: SNESLineSearchSet(snes,SNESLineSearchQuadratic,PETSC_NULL);
1138: break;
1139: case 3:
1140: SNESLineSearchSet(snes,SNESLineSearchCubic,PETSC_NULL);
1141: break;
1142: }
1143: }
1144: PetscOptionsTail();
1145: return(0);
1146: }
1147: /* -------------------------------------------------------------------------- */
1148: /*MC
1149: SNESLS - Newton based nonlinear solver that uses a line search
1151: Options Database:
1152: + -snes_ls [cubic,quadratic,basic,basicnonorms] - Selects line search
1153: . -snes_ls_alpha <alpha> - Sets alpha
1154: . -snes_ls_maxstep <maxstep> - Sets the maximum stepsize the line search will use (if the 2-norm(y) > maxstep then scale y to be y = (maxstep/2-norm(y)) *y)
1155: - -snes_ls_minlambda <minlambda> - Sets the minimum lambda the line search will use minlambda / max_i ( y[i]/x[i] )
1157: Notes: This is the default nonlinear solver in SNES
1159: Level: beginner
1161: .seealso: SNESCreate(), SNES, SNESSetType(), SNESTR, SNESLineSearchSet(),
1162: SNESLineSearchSetPostCheck(), SNESLineSearchNo(), SNESLineSearchCubic(), SNESLineSearchQuadratic(),
1163: SNESLineSearchSet(), SNESLineSearchNoNorms(), SNESLineSearchSetPreCheck()
1165: M*/
1169: PetscErrorCode SNESCreate_LS(SNES snes)
1170: {
1172: SNES_LS *neP;
1175: snes->ops->setup = SNESSetUp_LS;
1176: snes->ops->solve = SNESSolve_LS;
1177: snes->ops->destroy = SNESDestroy_LS;
1178: snes->ops->setfromoptions = SNESSetFromOptions_LS;
1179: snes->ops->view = SNESView_LS;
1181: PetscNewLog(snes,SNES_LS,&neP);
1182: snes->data = (void*)neP;
1183: neP->alpha = 1.e-4;
1184: neP->maxstep = 1.e8;
1185: neP->minlambda = 1.e-12;
1186: neP->LineSearch = SNESLineSearchCubic;
1187: neP->lsP = PETSC_NULL;
1188: neP->postcheckstep = PETSC_NULL;
1189: neP->postcheck = PETSC_NULL;
1190: neP->precheckstep = PETSC_NULL;
1191: neP->precheck = PETSC_NULL;
1193: PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSet_C",
1194: "SNESLineSearchSet_LS",
1195: SNESLineSearchSet_LS);
1196: PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSetPostCheck_C",
1197: "SNESLineSearchSetPostCheck_LS",
1198: SNESLineSearchSetPostCheck_LS);
1199: PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSetPreCheck_C",
1200: "SNESLineSearchSetPreCheck_LS",
1201: SNESLineSearchSetPreCheck_LS);
1203: return(0);
1204: }