1 #include <petsc/private/taolinesearchimpl.h>
2 #include <../src/tao/linesearch/impls/armijo/armijo.h>
3 
4 #define REPLACE_FIFO 1
5 #define REPLACE_MRU  2
6 
7 #define REFERENCE_MAX  1
8 #define REFERENCE_AVE  2
9 #define REFERENCE_MEAN 3
10 
TaoLineSearchDestroy_Armijo(TaoLineSearch ls)11 static PetscErrorCode TaoLineSearchDestroy_Armijo(TaoLineSearch ls)
12 {
13   TaoLineSearch_ARMIJO *armP = (TaoLineSearch_ARMIJO *)ls->data;
14   PetscErrorCode       ierr;
15 
16   PetscFunctionBegin;
17   ierr = PetscFree(armP->memory);CHKERRQ(ierr);
18   ierr = VecDestroy(&armP->x);CHKERRQ(ierr);
19   ierr = VecDestroy(&armP->work);CHKERRQ(ierr);
20   ierr = PetscFree(ls->data);CHKERRQ(ierr);
21   PetscFunctionReturn(0);
22 }
23 
TaoLineSearchReset_Armijo(TaoLineSearch ls)24 static PetscErrorCode TaoLineSearchReset_Armijo(TaoLineSearch ls)
25 {
26   TaoLineSearch_ARMIJO *armP = (TaoLineSearch_ARMIJO *)ls->data;
27   PetscErrorCode       ierr;
28 
29   PetscFunctionBegin;
30   ierr = PetscFree(armP->memory);CHKERRQ(ierr);
31   armP->memorySetup = PETSC_FALSE;
32   PetscFunctionReturn(0);
33 }
34 
TaoLineSearchSetFromOptions_Armijo(PetscOptionItems * PetscOptionsObject,TaoLineSearch ls)35 static PetscErrorCode TaoLineSearchSetFromOptions_Armijo(PetscOptionItems *PetscOptionsObject,TaoLineSearch ls)
36 {
37   TaoLineSearch_ARMIJO *armP = (TaoLineSearch_ARMIJO *)ls->data;
38   PetscErrorCode       ierr;
39 
40   PetscFunctionBegin;
41   ierr = PetscOptionsHead(PetscOptionsObject,"Armijo linesearch options");CHKERRQ(ierr);
42   ierr = PetscOptionsReal("-tao_ls_armijo_alpha", "initial reference constant", "", armP->alpha, &armP->alpha,NULL);CHKERRQ(ierr);
43   ierr = PetscOptionsReal("-tao_ls_armijo_beta_inf", "decrease constant one", "", armP->beta_inf, &armP->beta_inf,NULL);CHKERRQ(ierr);
44   ierr = PetscOptionsReal("-tao_ls_armijo_beta", "decrease constant", "", armP->beta, &armP->beta,NULL);CHKERRQ(ierr);
45   ierr = PetscOptionsReal("-tao_ls_armijo_sigma", "acceptance constant", "", armP->sigma, &armP->sigma,NULL);CHKERRQ(ierr);
46   ierr = PetscOptionsInt("-tao_ls_armijo_memory_size", "number of historical elements", "", armP->memorySize, &armP->memorySize,NULL);CHKERRQ(ierr);
47   ierr = PetscOptionsInt("-tao_ls_armijo_reference_policy", "policy for updating reference value", "", armP->referencePolicy, &armP->referencePolicy,NULL);CHKERRQ(ierr);
48   ierr = PetscOptionsInt("-tao_ls_armijo_replacement_policy", "policy for updating memory", "", armP->replacementPolicy, &armP->replacementPolicy,NULL);CHKERRQ(ierr);
49   ierr = PetscOptionsBool("-tao_ls_armijo_nondescending","Use nondescending armijo algorithm","",armP->nondescending,&armP->nondescending,NULL);CHKERRQ(ierr);
50   ierr = PetscOptionsTail();CHKERRQ(ierr);
51   PetscFunctionReturn(0);
52 }
53 
TaoLineSearchView_Armijo(TaoLineSearch ls,PetscViewer pv)54 static PetscErrorCode TaoLineSearchView_Armijo(TaoLineSearch ls, PetscViewer pv)
55 {
56   TaoLineSearch_ARMIJO *armP = (TaoLineSearch_ARMIJO *)ls->data;
57   PetscBool            isascii;
58   PetscErrorCode       ierr;
59 
60   PetscFunctionBegin;
61   ierr = PetscObjectTypeCompare((PetscObject)pv, PETSCVIEWERASCII, &isascii);CHKERRQ(ierr);
62   if (isascii) {
63     ierr=PetscViewerASCIIPrintf(pv,"  Armijo linesearch",armP->alpha);CHKERRQ(ierr);
64     if (armP->nondescending) {
65       ierr = PetscViewerASCIIPrintf(pv, " (nondescending)");CHKERRQ(ierr);
66     }
67     if (ls->bounded) {
68       ierr = PetscViewerASCIIPrintf(pv," (projected)");CHKERRQ(ierr);
69     }
70     ierr=PetscViewerASCIIPrintf(pv,": alpha=%g beta=%g ",(double)armP->alpha,(double)armP->beta);CHKERRQ(ierr);
71     ierr=PetscViewerASCIIPrintf(pv,"sigma=%g ",(double)armP->sigma);CHKERRQ(ierr);
72     ierr=PetscViewerASCIIPrintf(pv,"memsize=%D\n",armP->memorySize);CHKERRQ(ierr);
73   }
74   PetscFunctionReturn(0);
75 }
76 
77 /* @ TaoApply_Armijo - This routine performs a linesearch. It
78    backtracks until the (nonmonotone) Armijo conditions are satisfied.
79 
80    Input Parameters:
81 +  tao - Tao context
82 .  X - current iterate (on output X contains new iterate, X + step*S)
83 .  S - search direction
84 .  f - merit function evaluated at X
85 .  G - gradient of merit function evaluated at X
86 .  W - work vector
87 -  step - initial estimate of step length
88 
89    Output parameters:
90 +  f - merit function evaluated at new iterate, X + step*S
91 .  G - gradient of merit function evaluated at new iterate, X + step*S
92 .  X - new iterate
93 -  step - final step length
94 
95 @ */
TaoLineSearchApply_Armijo(TaoLineSearch ls,Vec x,PetscReal * f,Vec g,Vec s)96 static PetscErrorCode TaoLineSearchApply_Armijo(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, Vec s)
97 {
98   TaoLineSearch_ARMIJO *armP = (TaoLineSearch_ARMIJO *)ls->data;
99   PetscErrorCode       ierr;
100   PetscInt             i,its=0;
101   PetscReal            fact, ref, gdx;
102   PetscInt             idx;
103   PetscBool            g_computed=PETSC_FALSE; /* to prevent extra gradient computation */
104 
105   PetscFunctionBegin;
106   ierr = TaoLineSearchMonitor(ls, 0, *f, 0.0);CHKERRQ(ierr);
107 
108   ls->reason = TAOLINESEARCH_CONTINUE_ITERATING;
109   if (!armP->work) {
110     ierr = VecDuplicate(x,&armP->work);CHKERRQ(ierr);
111     armP->x = x;
112     ierr = PetscObjectReference((PetscObject)armP->x);CHKERRQ(ierr);
113   } else if (x != armP->x) {
114     /* If x has changed, then recreate work */
115     ierr = VecDestroy(&armP->work);CHKERRQ(ierr);
116     ierr = VecDuplicate(x,&armP->work);CHKERRQ(ierr);
117     ierr = PetscObjectDereference((PetscObject)armP->x);CHKERRQ(ierr);
118     armP->x = x;
119     ierr = PetscObjectReference((PetscObject)armP->x);CHKERRQ(ierr);
120   }
121 
122   /* Check linesearch parameters */
123   if (armP->alpha < 1) {
124     ierr = PetscInfo1(ls,"Armijo line search error: alpha (%g) < 1\n", (double)armP->alpha);CHKERRQ(ierr);
125     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
126   } else if ((armP->beta <= 0) || (armP->beta >= 1)) {
127     ierr = PetscInfo1(ls,"Armijo line search error: beta (%g) invalid\n", (double)armP->beta);CHKERRQ(ierr);
128     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
129   } else if ((armP->beta_inf <= 0) || (armP->beta_inf >= 1)) {
130     ierr = PetscInfo1(ls,"Armijo line search error: beta_inf (%g) invalid\n", (double)armP->beta_inf);CHKERRQ(ierr);
131     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
132   } else if ((armP->sigma <= 0) || (armP->sigma >= 0.5)) {
133     ierr = PetscInfo1(ls,"Armijo line search error: sigma (%g) invalid\n", (double)armP->sigma);CHKERRQ(ierr);
134     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
135   } else if (armP->memorySize < 1) {
136     ierr = PetscInfo1(ls,"Armijo line search error: memory_size (%D) < 1\n", armP->memorySize);CHKERRQ(ierr);
137     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
138   } else if ((armP->referencePolicy != REFERENCE_MAX) && (armP->referencePolicy != REFERENCE_AVE) && (armP->referencePolicy != REFERENCE_MEAN)) {
139     ierr = PetscInfo(ls,"Armijo line search error: reference_policy invalid\n");CHKERRQ(ierr);
140     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
141   } else if ((armP->replacementPolicy != REPLACE_FIFO) && (armP->replacementPolicy != REPLACE_MRU)) {
142     ierr = PetscInfo(ls,"Armijo line search error: replacement_policy invalid\n");CHKERRQ(ierr);
143     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
144   } else if (PetscIsInfOrNanReal(*f)) {
145     ierr = PetscInfo(ls,"Armijo line search error: initial function inf or nan\n");CHKERRQ(ierr);
146     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
147   }
148 
149   if (ls->reason != TAOLINESEARCH_CONTINUE_ITERATING) {
150     PetscFunctionReturn(0);
151   }
152 
153   /* Check to see of the memory has been allocated.  If not, allocate
154      the historical array and populate it with the initial function
155      values. */
156   if (!armP->memory) {
157     ierr = PetscMalloc1(armP->memorySize, &armP->memory);CHKERRQ(ierr);
158   }
159 
160   if (!armP->memorySetup) {
161     for (i = 0; i < armP->memorySize; i++) {
162       armP->memory[i] = armP->alpha*(*f);
163     }
164 
165     armP->current = 0;
166     armP->lastReference = armP->memory[0];
167     armP->memorySetup=PETSC_TRUE;
168   }
169 
170   /* Calculate reference value (MAX) */
171   ref = armP->memory[0];
172   idx = 0;
173 
174   for (i = 1; i < armP->memorySize; i++) {
175     if (armP->memory[i] > ref) {
176       ref = armP->memory[i];
177       idx = i;
178     }
179   }
180 
181   if (armP->referencePolicy == REFERENCE_AVE) {
182     ref = 0;
183     for (i = 0; i < armP->memorySize; i++) {
184       ref += armP->memory[i];
185     }
186     ref = ref / armP->memorySize;
187     ref = PetscMax(ref, armP->memory[armP->current]);
188   } else if (armP->referencePolicy == REFERENCE_MEAN) {
189     ref = PetscMin(ref, 0.5*(armP->lastReference + armP->memory[armP->current]));
190   }
191   ierr = VecDot(g,s,&gdx);CHKERRQ(ierr);
192 
193   if (PetscIsInfOrNanReal(gdx)) {
194     ierr = PetscInfo1(ls,"Initial Line Search step * g is Inf or Nan (%g)\n",(double)gdx);CHKERRQ(ierr);
195     ls->reason=TAOLINESEARCH_FAILED_INFORNAN;
196     PetscFunctionReturn(0);
197   }
198   if (gdx >= 0.0) {
199     ierr = PetscInfo1(ls,"Initial Line Search step is not descent direction (g's=%g)\n",(double)gdx);CHKERRQ(ierr);
200     ls->reason = TAOLINESEARCH_FAILED_ASCENT;
201     PetscFunctionReturn(0);
202   }
203 
204   if (armP->nondescending) {
205     fact = armP->sigma;
206   } else {
207     fact = armP->sigma * gdx;
208   }
209   ls->step = ls->initstep;
210   while (ls->step >= ls->stepmin && (ls->nfeval+ls->nfgeval) < ls->max_funcs) {
211     /* Calculate iterate */
212     ++its;
213     ierr = VecCopy(x,armP->work);CHKERRQ(ierr);
214     ierr = VecAXPY(armP->work,ls->step,s);CHKERRQ(ierr);
215     if (ls->bounded) {
216       ierr = VecMedian(ls->lower,armP->work,ls->upper,armP->work);CHKERRQ(ierr);
217     }
218 
219     /* Calculate function at new iterate */
220     if (ls->hasobjective) {
221       ierr = TaoLineSearchComputeObjective(ls,armP->work,f);CHKERRQ(ierr);
222       g_computed=PETSC_FALSE;
223     } else if (ls->usegts) {
224       ierr = TaoLineSearchComputeObjectiveAndGTS(ls,armP->work,f,&gdx);CHKERRQ(ierr);
225       g_computed=PETSC_FALSE;
226     } else {
227       ierr = TaoLineSearchComputeObjectiveAndGradient(ls,armP->work,f,g);CHKERRQ(ierr);
228       g_computed=PETSC_TRUE;
229     }
230     if (ls->step == ls->initstep) {
231       ls->f_fullstep = *f;
232     }
233 
234     ierr = TaoLineSearchMonitor(ls, its, *f, ls->step);CHKERRQ(ierr);
235 
236     if (PetscIsInfOrNanReal(*f)) {
237       ls->step *= armP->beta_inf;
238     } else {
239       /* Check descent condition */
240       if (armP->nondescending && *f <= ref - ls->step*fact*ref)
241         break;
242       if (!armP->nondescending && *f <= ref + ls->step*fact) {
243         break;
244       }
245 
246       ls->step *= armP->beta;
247     }
248   }
249 
250   /* Check termination */
251   if (PetscIsInfOrNanReal(*f)) {
252     ierr = PetscInfo(ls, "Function is inf or nan.\n");CHKERRQ(ierr);
253     ls->reason = TAOLINESEARCH_FAILED_INFORNAN;
254   } else if (ls->step < ls->stepmin) {
255     ierr = PetscInfo(ls, "Step length is below tolerance.\n");CHKERRQ(ierr);
256     ls->reason = TAOLINESEARCH_HALTED_RTOL;
257   } else if ((ls->nfeval+ls->nfgeval) >= ls->max_funcs) {
258     ierr = PetscInfo2(ls, "Number of line search function evals (%D) > maximum allowed (%D)\n",ls->nfeval+ls->nfgeval, ls->max_funcs);CHKERRQ(ierr);
259     ls->reason = TAOLINESEARCH_HALTED_MAXFCN;
260   }
261   if (ls->reason) {
262     PetscFunctionReturn(0);
263   }
264 
265   /* Successful termination, update memory */
266   ls->reason = TAOLINESEARCH_SUCCESS;
267   armP->lastReference = ref;
268   if (armP->replacementPolicy == REPLACE_FIFO) {
269     armP->memory[armP->current++] = *f;
270     if (armP->current >= armP->memorySize) {
271       armP->current = 0;
272     }
273   } else {
274     armP->current = idx;
275     armP->memory[idx] = *f;
276   }
277 
278   /* Update iterate and compute gradient */
279   ierr = VecCopy(armP->work,x);CHKERRQ(ierr);
280   if (!g_computed) {
281     ierr = TaoLineSearchComputeGradient(ls, x, g);CHKERRQ(ierr);
282   }
283   ierr = PetscInfo2(ls, "%D function evals in line search, step = %g\n",ls->nfeval, (double)ls->step);CHKERRQ(ierr);
284   PetscFunctionReturn(0);
285 }
286 
287 /*MC
288    TAOLINESEARCHARMIJO - Backtracking line-search that satisfies only the (nonmonotone) Armijo condition
289    (i.e., sufficient decrease).
290 
291    Armijo line-search type can be selected with "-tao_ls_type armijo".
292 
293    Level: developer
294 
295 .seealso: TaoLineSearchCreate(), TaoLineSearchSetType(), TaoLineSearchApply()
296 
297 .keywords: Tao, linesearch
298 M*/
TaoLineSearchCreate_Armijo(TaoLineSearch ls)299 PETSC_EXTERN PetscErrorCode TaoLineSearchCreate_Armijo(TaoLineSearch ls)
300 {
301   TaoLineSearch_ARMIJO *armP;
302   PetscErrorCode       ierr;
303 
304   PetscFunctionBegin;
305   PetscValidHeaderSpecific(ls,TAOLINESEARCH_CLASSID,1);
306   ierr = PetscNewLog(ls,&armP);CHKERRQ(ierr);
307 
308   armP->memory = NULL;
309   armP->alpha = 1.0;
310   armP->beta = 0.5;
311   armP->beta_inf = 0.5;
312   armP->sigma = 1e-4;
313   armP->memorySize = 1;
314   armP->referencePolicy = REFERENCE_MAX;
315   armP->replacementPolicy = REPLACE_MRU;
316   armP->nondescending=PETSC_FALSE;
317   ls->data = (void*)armP;
318   ls->initstep=1.0;
319   ls->ops->setup = NULL;
320   ls->ops->apply = TaoLineSearchApply_Armijo;
321   ls->ops->view = TaoLineSearchView_Armijo;
322   ls->ops->destroy = TaoLineSearchDestroy_Armijo;
323   ls->ops->reset = TaoLineSearchReset_Armijo;
324   ls->ops->setfromoptions = TaoLineSearchSetFromOptions_Armijo;
325   PetscFunctionReturn(0);
326 }
327