xref: /petsc/src/snes/linesearch/impls/cp/linesearchcp.c (revision 1db2f0e57964f94f63a3f5cf6bb1e5befdb56f30)
1 #include <private/linesearchimpl.h>
2 #include <petscsnes.h>
3 
4 #undef __FUNCT__
5 #define __FUNCT__ "LineSearchApply_CP"
6 
7 /*@C
8    LineSearchCP - This routine is not a line search at all;
9    it simply uses the full step.  Thus, this routine is intended
10    to serve as a template and is not recommended for general use.
11 
12    Logically Collective on SNES and Vec
13 
14    Input Parameters:
15 +  snes - nonlinear context
16 .  lsctx - optional context for line search (not used here)
17 .  x - current iterate
18 .  f - residual evaluated at x
19 .  y - search direction
20 .  fnorm - 2-norm of f
21 -  xnorm - norm of x if known, otherwise 0
22 
23    Output Parameters:
24 +  g - residual evaluated at new iterate y
25 .  w - new iterate
26 .  gnorm - 2-norm of g
27 .  ynorm - 2-norm of search length
28 -  flag - PETSC_TRUE on success, PETSC_FALSE on failure
29 
30    Options Database Key:
31 .  -snes_ls basic - Activates SNESLineSearchNo()
32 
33    Level: advanced
34 
35 .keywords: SNES, nonlinear, line search, cubic
36 
37 .seealso: SNESLineSearchCubic(), SNESLineSearchQuadratic(),
38           SNESLineSearchSet(), SNESLineSearchNoNorms()
39 @*/
40 PetscErrorCode  LineSearchApply_CP(LineSearch linesearch)
41 {
42 
43   PetscBool      changed_y, changed_w, domainerror;
44   PetscErrorCode ierr;
45   Vec             X, Y, F, W;
46   SNES            snes;
47   PetscReal       xnorm, ynorm, gnorm, steptol, maxstep;
48 
49   PetscReal       lambda, lambda_old, lambda_update, delLambda;
50   PetscScalar     fty, fty_old;
51   PetscInt        i, max_its;
52 
53   PetscViewer     monitor;
54 
55   PetscFunctionBegin;
56 
57   ierr = LineSearchGetVecs(linesearch, &X, &F, &Y, &W, PETSC_NULL);CHKERRQ(ierr);
58   ierr = LineSearchGetNorms(linesearch, &xnorm, &gnorm, &ynorm);CHKERRQ(ierr);
59   ierr = LineSearchGetSNES(linesearch, &snes);CHKERRQ(ierr);
60   ierr = LineSearchGetLambda(linesearch, &lambda);CHKERRQ(ierr);
61   ierr = LineSearchGetMaxIts(linesearch, &max_its);CHKERRQ(ierr);
62   ierr = LineSearchGetStepTolerance(linesearch, &steptol);CHKERRQ(ierr);
63   ierr = LineSearchGetMaxStep(linesearch, &maxstep);CHKERRQ(ierr);
64   ierr = LineSearchSetSuccess(linesearch, PETSC_TRUE);CHKERRQ(ierr);
65   ierr = LineSearchGetMonitor(linesearch, &monitor);CHKERRQ(ierr);
66 
67   /* precheck */
68   ierr = LineSearchPreCheck(linesearch, &changed_y);CHKERRQ(ierr);
69   lambda_old = 0.0;
70   ierr = VecDot(F, Y, &fty_old);CHKERRQ(ierr);
71   for (i = 0; i < max_its; i++) {
72 
73     /* compute the norm at lambda */
74     ierr = VecCopy(X, W);CHKERRQ(ierr);
75     ierr = VecAXPY(W, -lambda, Y);CHKERRQ(ierr);
76     ierr = SNESComputeFunction(snes, W, F);CHKERRQ(ierr);
77 
78     ierr = VecDot(F, Y, &fty);CHKERRQ(ierr);
79 
80     delLambda    = lambda - lambda_old;
81     if (PetscAbsReal(delLambda) < steptol) break;
82     if (monitor) {
83       ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
84       ierr = PetscViewerASCIIPrintf(monitor,"    Line search: lambdas = [%g, %g], ftys = [%g, %g]\n",
85                                     lambda, lambda_old, PetscRealPart(fty), PetscRealPart(fty_old));CHKERRQ(ierr);
86       ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
87     }
88 
89     /* compute the search direction */
90     lambda_update =  PetscRealPart((fty*lambda_old - fty_old*lambda) / (fty - fty_old));
91     if (PetscIsInfOrNanScalar(lambda_update)) break;
92     if (lambda_update > maxstep) {
93       break;
94     }
95 
96     /* compute the new state of the line search */
97     lambda_old = lambda;
98     lambda = lambda_update;
99     fty_old = fty;
100   }
101   /* construct the solution */
102   ierr = VecCopy(X, W);CHKERRQ(ierr);
103   ierr = VecAXPY(W, -lambda, Y);CHKERRQ(ierr);
104 
105   /* postcheck */
106   ierr = LineSearchPostCheck(linesearch, &changed_y, &changed_w);CHKERRQ(ierr);
107   if (changed_y) {
108     ierr = VecAXPY(X, -lambda, Y);CHKERRQ(ierr);
109   } else {
110     ierr = VecCopy(W, X);CHKERRQ(ierr);
111   }
112   ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr);
113   ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr);
114   if (domainerror) {
115     ierr = LineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
116     PetscFunctionReturn(0);
117   }
118 
119   ierr = LineSearchComputeNorms(linesearch);CHKERRQ(ierr);
120   ierr = LineSearchGetNorms(linesearch, &xnorm, &gnorm, &ynorm);CHKERRQ(ierr);
121 
122   ierr = LineSearchSetLambda(linesearch, lambda);CHKERRQ(ierr);
123 
124   if (monitor) {
125     ierr = PetscViewerASCIIAddTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
126     ierr = PetscViewerASCIIPrintf(monitor,"    Line search terminated: lambda = %g, fnorms = %g\n", lambda, gnorm);CHKERRQ(ierr);
127     ierr = PetscViewerASCIISubtractTab(monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
128   }
129   if (lambda <= steptol) {
130     ierr = LineSearchSetSuccess(linesearch, PETSC_FALSE);CHKERRQ(ierr);
131   }
132   PetscFunctionReturn(0);
133 }
134 
135 EXTERN_C_BEGIN
136 #undef __FUNCT__
137 #define __FUNCT__ "LineSearchCreate_CP"
138 PetscErrorCode LineSearchCreate_CP(LineSearch linesearch)
139 {
140   PetscFunctionBegin;
141   linesearch->ops->apply          = LineSearchApply_CP;
142   linesearch->ops->destroy        = PETSC_NULL;
143   linesearch->ops->setfromoptions = PETSC_NULL;
144   linesearch->ops->reset          = PETSC_NULL;
145   linesearch->ops->view           = PETSC_NULL;
146   linesearch->ops->setup          = PETSC_NULL;
147   PetscFunctionReturn(0);
148 }
149 EXTERN_C_END
150