1 2 #include <petsc/private/taolinesearchimpl.h> 3 4 static PetscErrorCode TaoLineSearchDestroy_Unit(TaoLineSearch ls) 5 { 6 PetscFunctionBegin; 7 PetscCall(PetscFree(ls->data)); 8 PetscFunctionReturn(0); 9 } 10 11 static PetscErrorCode TaoLineSearchSetFromOptions_Unit(PetscOptionItems *PetscOptionsObject,TaoLineSearch ls) 12 { 13 PetscFunctionBegin; 14 PetscCall(PetscOptionsHead(PetscOptionsObject,"No Unit line search options")); 15 PetscCall(PetscOptionsTail()); 16 PetscFunctionReturn(0); 17 } 18 19 static PetscErrorCode TaoLineSearchView_Unit(TaoLineSearch ls,PetscViewer viewer) 20 { 21 PetscBool isascii; 22 23 PetscFunctionBegin; 24 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii)); 25 if (isascii) { 26 PetscCall(PetscViewerASCIIPrintf(viewer," Line Search: Unit Step.\n")); 27 } 28 PetscFunctionReturn(0); 29 } 30 31 static PetscErrorCode TaoLineSearchApply_Unit(TaoLineSearch ls,Vec x,PetscReal *f,Vec g,Vec step_direction) 32 { 33 PetscReal ftry; 34 PetscReal startf = *f; 35 36 PetscFunctionBegin; 37 /* Take unit step (newx = startx + 1.0*step_direction) */ 38 PetscCall(TaoLineSearchMonitor(ls, 0, *f, 0.0)); 39 PetscCall(VecAXPY(x,1.0,step_direction)); 40 PetscCall(TaoLineSearchComputeObjectiveAndGradient(ls,x,&ftry,g)); 41 PetscCall(TaoLineSearchMonitor(ls, 1, *f, 1.0)); 42 PetscCall(PetscInfo(ls,"Tao Apply Unit Step: %4.4e\n",1.0)); 43 if (startf < ftry) { 44 PetscCall(PetscInfo(ls,"Tao Apply Unit Step, FINCREASE: F old:= %12.10e, F new: %12.10e\n",(double)startf,(double)ftry)); 45 } 46 *f = ftry; 47 ls->step = 1.0; 48 ls->reason=TAOLINESEARCH_SUCCESS; 49 PetscFunctionReturn(0); 50 } 51 52 /*MC 53 TAOLINESEARCHUNIT - Line-search type that disables line search and accepts the unit step length every time 54 55 Level: developer 56 57 .seealso: TaoLineSearchCreate(), TaoLineSearchSetType(), TaoLineSearchApply() 58 59 .keywords: Tao, linesearch 60 M*/ 61 PETSC_EXTERN PetscErrorCode TaoLineSearchCreate_Unit(TaoLineSearch ls) 62 { 63 PetscFunctionBegin; 64 ls->ops->setup = NULL; 65 ls->ops->reset = NULL; 66 ls->ops->monitor = NULL; 67 ls->ops->apply = TaoLineSearchApply_Unit; 68 ls->ops->view = TaoLineSearchView_Unit; 69 ls->ops->destroy = TaoLineSearchDestroy_Unit; 70 ls->ops->setfromoptions = TaoLineSearchSetFromOptions_Unit; 71 PetscFunctionReturn(0); 72 } 73