1 #include <private/linesearchimpl.h> 2 #include <private/snesimpl.h> 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "LineSearchApply_Basic" 6 7 /*@C 8 SNESLineSearchBasic - 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_Basic(LineSearch linesearch) 41 { 42 PetscBool changed_y, changed_w; 43 PetscErrorCode ierr; 44 Vec X, F, Y, W; 45 SNES snes; 46 PetscReal gnorm, xnorm, ynorm, lambda; 47 PetscBool domainerror; 48 49 PetscFunctionBegin; 50 51 ierr = LineSearchGetVecs(linesearch, &X, &F, &Y, &W, PETSC_NULL);CHKERRQ(ierr); 52 ierr = LineSearchGetNorms(linesearch, &xnorm, &gnorm, &ynorm);CHKERRQ(ierr); 53 ierr = LineSearchGetLambda(linesearch, &lambda);CHKERRQ(ierr); 54 ierr = LineSearchGetSNES(linesearch, &snes);CHKERRQ(ierr); 55 ierr = LineSearchSetSuccess(linesearch, PETSC_TRUE);CHKERRQ(ierr); 56 57 /* precheck */ 58 ierr = LineSearchPreCheck(linesearch, &changed_y);CHKERRQ(ierr); 59 60 /* update */ 61 ierr = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr); 62 63 /* postcheck */ 64 ierr = LineSearchPostCheck(linesearch, &changed_y, &changed_w);CHKERRQ(ierr); 65 if (changed_y) { 66 ierr = VecWAXPY(W,-lambda,Y,X);CHKERRQ(ierr); 67 } 68 ierr = SNESComputeFunction(snes,W,F);CHKERRQ(ierr); 69 ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr); 70 if (domainerror) { 71 ierr = LineSearchSetSuccess(linesearch, PETSC_FALSE); 72 PetscFunctionReturn(0); 73 } 74 75 ierr = LineSearchComputeNorms(linesearch);CHKERRQ(ierr); 76 77 /* copy the solution over */ 78 ierr = VecCopy(W, X);CHKERRQ(ierr); 79 PetscFunctionReturn(0); 80 } 81 82 EXTERN_C_BEGIN 83 #undef __FUNCT__ 84 #define __FUNCT__ "LineSearchCreate_Basic" 85 PetscErrorCode LineSearchCreate_Basic(LineSearch linesearch) 86 { 87 PetscFunctionBegin; 88 linesearch->ops->apply = LineSearchApply_Basic; 89 linesearch->ops->destroy = PETSC_NULL; 90 linesearch->ops->setfromoptions = PETSC_NULL; 91 linesearch->ops->reset = PETSC_NULL; 92 linesearch->ops->view = PETSC_NULL; 93 linesearch->ops->setup = PETSC_NULL; 94 PetscFunctionReturn(0); 95 } 96 EXTERN_C_END 97