1 #include <petsc/private/taolinesearchimpl.h> 2 3 static PetscErrorCode TaoLineSearchDestroy_Unit(TaoLineSearch ls) 4 { 5 PetscFunctionBegin; 6 PetscFunctionReturn(PETSC_SUCCESS); 7 } 8 9 static PetscErrorCode TaoLineSearchSetFromOptions_Unit(TaoLineSearch ls, PetscOptionItems *PetscOptionsObject) 10 { 11 PetscFunctionBegin; 12 PetscFunctionReturn(PETSC_SUCCESS); 13 } 14 15 static PetscErrorCode TaoLineSearchView_Unit(TaoLineSearch ls, PetscViewer viewer) 16 { 17 PetscBool isascii; 18 19 PetscFunctionBegin; 20 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii)); 21 if (isascii) PetscCall(PetscViewerASCIIPrintf(viewer, " Line Search: Unit Step %g.\n", (double)ls->initstep)); 22 PetscFunctionReturn(PETSC_SUCCESS); 23 } 24 25 /* Take unit step (newx = startx + initstep*step_direction) */ 26 static PetscErrorCode TaoLineSearchApply_Unit(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, Vec step_direction) 27 { 28 PetscFunctionBegin; 29 PetscCall(TaoLineSearchMonitor(ls, 0, *f, 0.0)); 30 ls->step = ls->initstep; 31 PetscCall(VecAXPY(x, ls->step, step_direction)); 32 PetscCall(TaoLineSearchComputeObjectiveAndGradient(ls, x, f, g)); 33 PetscCall(TaoLineSearchMonitor(ls, 1, *f, ls->step)); 34 ls->reason = TAOLINESEARCH_SUCCESS; 35 PetscFunctionReturn(PETSC_SUCCESS); 36 } 37 38 /*MC 39 TAOLINESEARCHUNIT - Line-search type that disables line search and accepts the unit step length every time 40 41 Options Database Keys: 42 . -tao_ls_stepinit <step> - steplength 43 44 Level: developer 45 46 .seealso: `Tao`, `TaoLinesearch`, `TaoLineSearchCreate()`, `TaoLineSearchSetType()`, `TaoLineSearchApply()` 47 M*/ 48 PETSC_EXTERN PetscErrorCode TaoLineSearchCreate_Unit(TaoLineSearch ls) 49 { 50 PetscFunctionBegin; 51 ls->ops->setup = NULL; 52 ls->ops->reset = NULL; 53 ls->ops->monitor = NULL; 54 ls->ops->apply = TaoLineSearchApply_Unit; 55 ls->ops->view = TaoLineSearchView_Unit; 56 ls->ops->destroy = TaoLineSearchDestroy_Unit; 57 ls->ops->setfromoptions = TaoLineSearchSetFromOptions_Unit; 58 PetscFunctionReturn(PETSC_SUCCESS); 59 } 60