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