xref: /petsc/src/snes/linesearch/impls/cp/linesearchcp.c (revision bbd5d0b317af04eb83c39e72581bf86c64bd741a)
1 #include <private/linesearchimpl.h>
2 #include <private/snesimpl.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;
44   PetscErrorCode ierr;
45   Vec             X = linesearch->vec_sol;
46   Vec             F  = linesearch->vec_func;
47   Vec             Y  = linesearch->vec_update;
48   Vec             W =  linesearch->vec_sol_new;
49   SNES            snes  = linesearch->snes;
50   PetscReal       *gnorm  = &linesearch->fnorm;
51   PetscReal       *ynorm  = &linesearch->ynorm;
52   PetscReal       *xnorm =  &linesearch->xnorm;
53 
54   PetscReal       lambda, lambda_old, lambda_update, delLambda;
55   PetscScalar     fty, fty_old;
56   PetscInt        i;
57 
58   PetscFunctionBegin;
59   /* precheck */
60   ierr = LineSearchPreCheck(linesearch, &changed_y);CHKERRQ(ierr);
61   lambda = linesearch->lambda;
62   lambda_old = 0.0;
63   ierr = VecDot(F, Y, &fty_old);CHKERRQ(ierr);
64   linesearch->success = PETSC_TRUE;
65   for (i = 0; i < linesearch->max_its; i++) {
66 
67     /* compute the norm at lambda */
68     ierr = VecCopy(X, W);CHKERRQ(ierr);
69     ierr = VecAXPY(W, -lambda, Y);CHKERRQ(ierr);
70     ierr = SNESComputeFunction(snes, W, F);CHKERRQ(ierr);
71 
72     ierr = VecDot(F, Y, &fty);CHKERRQ(ierr);
73 
74     delLambda    = lambda - lambda_old;
75     if (PetscAbsReal(delLambda) < linesearch->steptol) break;
76     if (linesearch->monitor) {
77       ierr = PetscViewerASCIIAddTab(linesearch->monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
78       ierr = PetscViewerASCIIPrintf(linesearch->monitor,"    Line search: lambdas = [%g, %g], ftys = [%g, %g]\n",
79                                     lambda, lambda_old, PetscRealPart(fty), PetscRealPart(fty_old));CHKERRQ(ierr);
80       ierr = PetscViewerASCIISubtractTab(linesearch->monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
81     }
82 
83     /* compute the search direction */
84     lambda_update =  PetscRealPart((fty*lambda_old - fty_old*lambda) / (fty - fty_old));
85     if (PetscIsInfOrNanScalar(lambda_update)) break;
86     if (lambda_update > linesearch->maxstep) {
87       break;
88     }
89 
90     /* compute the new state of the line search */
91     lambda_old = lambda;
92     lambda = lambda_update;
93     fty_old = fty;
94   }
95   /* construct the solution */
96   ierr = VecCopy(X, W);CHKERRQ(ierr);
97   ierr = VecAXPY(W, -lambda, Y);CHKERRQ(ierr);
98 
99   /* postcheck */
100   ierr = LineSearchPostCheck(linesearch, &changed_y, &changed_w);CHKERRQ(ierr);
101   if (changed_y) {
102     ierr = VecAXPY(X, -lambda, Y);CHKERRQ(ierr);
103   } else {
104     ierr = VecCopy(W, X);CHKERRQ(ierr);
105   }
106   ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr);
107   if (snes->domainerror) {
108     linesearch->success = PETSC_FALSE;
109     PetscFunctionReturn(0);
110   }
111 
112   ierr = VecNormBegin(F, NORM_2, gnorm);CHKERRQ(ierr);
113   ierr = VecNormBegin(X, NORM_2, xnorm);CHKERRQ(ierr);
114   ierr = VecNormBegin(Y, NORM_2, ynorm);CHKERRQ(ierr);
115   ierr = VecNormEnd(F, NORM_2, gnorm);CHKERRQ(ierr);
116   ierr = VecNormEnd(X, NORM_2, xnorm);CHKERRQ(ierr);
117   ierr = VecNormEnd(Y, NORM_2, ynorm);CHKERRQ(ierr);
118 
119   linesearch->lambda = lambda;
120   if (linesearch->monitor) {
121     ierr = PetscViewerASCIIAddTab(linesearch->monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
122     ierr = PetscViewerASCIIPrintf(linesearch->monitor,"    Line search terminated: lambda = %g, fnorms = %g\n", lambda, linesearch->fnorm);CHKERRQ(ierr);
123     ierr = PetscViewerASCIISubtractTab(linesearch->monitor,((PetscObject)linesearch)->tablevel);CHKERRQ(ierr);
124   }
125   if (lambda <= linesearch->steptol) {
126     linesearch->success = PETSC_FALSE;
127   }
128   PetscFunctionReturn(0);
129 }
130 
131 EXTERN_C_BEGIN
132 #undef __FUNCT__
133 #define __FUNCT__ "LineSearchCreate_CP"
134 PetscErrorCode LineSearchCreate_CP(LineSearch linesearch)
135 {
136   PetscFunctionBegin;
137   linesearch->ops->apply          = LineSearchApply_CP;
138   linesearch->ops->destroy        = PETSC_NULL;
139   linesearch->ops->setfromoptions = PETSC_NULL;
140   linesearch->ops->reset          = PETSC_NULL;
141   linesearch->ops->view           = PETSC_NULL;
142   linesearch->ops->setup          = PETSC_NULL;
143   PetscFunctionReturn(0);
144 }
145 EXTERN_C_END
146