1 #include <petsc/private/linesearchimpl.h>
2 #include <petsc/private/snesimpl.h>
3
SNESLineSearchApply_Basic(SNESLineSearch linesearch)4 static PetscErrorCode SNESLineSearchApply_Basic(SNESLineSearch linesearch)
5 {
6 PetscBool changed_y, changed_w;
7 Vec X, F, Y, W;
8 SNES snes;
9 PetscReal gnorm, xnorm, ynorm, lambda, fnorm = 0.0;
10
11 PetscFunctionBegin;
12 PetscCall(SNESLineSearchGetVecs(linesearch, &X, &F, &Y, &W, NULL));
13 PetscCall(SNESLineSearchGetNorms(linesearch, &xnorm, &gnorm, &ynorm));
14 PetscCall(SNESLineSearchGetLambda(linesearch, &lambda));
15 PetscCall(SNESLineSearchGetSNES(linesearch, &snes));
16 PetscCall(SNESLineSearchSetReason(linesearch, SNES_LINESEARCH_SUCCEEDED));
17
18 /* precheck */
19 PetscCall(SNESLineSearchPreCheck(linesearch, X, Y, &changed_y));
20
21 /* update */
22 PetscCall(VecWAXPY(W, -lambda, Y, X));
23 if (linesearch->ops->viproject) PetscCall((*linesearch->ops->viproject)(snes, W));
24
25 /* postcheck */
26 PetscCall(SNESLineSearchPostCheck(linesearch, X, Y, W, &changed_y, &changed_w));
27 if (changed_y) {
28 if (!changed_w) PetscCall(VecWAXPY(W, -lambda, Y, X));
29 if (linesearch->ops->viproject) PetscCall((*linesearch->ops->viproject)(snes, W));
30 }
31 if (linesearch->norms || snes->iter < snes->max_its - 1) {
32 PetscCall((*linesearch->ops->snesfunc)(snes, W, F));
33 PetscCall(VecNorm(F, NORM_2, &fnorm));
34 SNESLineSearchCheckFunctionDomainError(snes, linesearch, fnorm);
35 }
36 if (linesearch->norms) {
37 PetscCall(VecNormBegin(Y, NORM_2, &linesearch->ynorm));
38 PetscCall(VecNormBegin(W, NORM_2, &linesearch->xnorm));
39 PetscCall(VecNormEnd(Y, NORM_2, &linesearch->ynorm));
40 PetscCall(VecNormEnd(W, NORM_2, &linesearch->xnorm));
41
42 if (linesearch->ops->vinorm) {
43 linesearch->fnorm = gnorm;
44
45 PetscCall((*linesearch->ops->vinorm)(snes, F, W, &linesearch->fnorm));
46 } else linesearch->fnorm = fnorm;
47 }
48
49 /* copy the solution over */
50 PetscCall(VecCopy(W, X));
51 PetscFunctionReturn(PETSC_SUCCESS);
52 }
53
54 /*MC
55 SNESLINESEARCHBASIC - This line search implementation is not a line
56 search at all; it simply uses the full step $x_{k+1} = x_k - \lambda Y_k$ with $\lambda=1$.
57 Alternatively, $\lambda$ can be configured to be a constant damping factor by setting `snes_linesearch_damping`.
58 Thus, this routine is intended for methods with well-scaled updates; i.e. Newton's method (`SNESNEWTONLS`), on
59 well-behaved problems. Also named `SNESLINESEARCHNONE`.
60
61 Options Database Keys:
62 + -snes_linesearch_damping <1.0> - step length is scaled by this factor
63 - -snes_linesearch_norms <true> - whether to compute norms or not (`SNESLineSearchSetComputeNorms()`)
64
65 Note:
66 For methods with ill-scaled updates (`SNESNRICHARDSON`, `SNESNCG`), a small
67 damping parameter may yield satisfactory, but slow convergence, despite
68 the lack of the line search.
69
70 Level: advanced
71
72 .seealso: [](ch_snes), `SNES`, `SNESLineSearch`, `SNESLineSearchType`, `SNESGetLineSearch()`, `SNESLineSearchCreate()`, `SNESLineSearchSetType()`, `SNESLineSearchSetDamping()`, `SNESLineSearchSetComputeNorms()`
73 M*/
SNESLineSearchCreate_Basic(SNESLineSearch linesearch)74 PETSC_EXTERN PetscErrorCode SNESLineSearchCreate_Basic(SNESLineSearch linesearch)
75 {
76 PetscFunctionBegin;
77 linesearch->ops->apply = SNESLineSearchApply_Basic;
78 linesearch->ops->destroy = NULL;
79 linesearch->ops->setfromoptions = NULL;
80 linesearch->ops->reset = NULL;
81 linesearch->ops->view = NULL;
82 linesearch->ops->setup = NULL;
83 PetscFunctionReturn(PETSC_SUCCESS);
84 }
85