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