xref: /petsc/src/tao/linesearch/impls/owarmijo/owarmijo.c (revision 87f595a510d358797c1cc8101e6f6c6930cb072f)
1 #include "petscvec.h"
2 #include "taosolver.h"
3 #include "tao-private/taolinesearch_impl.h"
4 #include "owarmijo.h"
5 
6 #define REPLACE_FIFO 1
7 #define REPLACE_MRU  2
8 
9 #define REFERENCE_MAX  1
10 #define REFERENCE_AVE  2
11 #define REFERENCE_MEAN 3
12 
13 
14 
15 static PetscErrorCode ProjWork_OWLQN(Vec w,Vec x,Vec gv,PetscReal *gdx)
16 {
17   PetscReal *xptr,*wptr,*gptr;
18   PetscErrorCode ierr;
19   PetscInt low,high,low1,high1,low2,high2,i;
20 
21   PetscFunctionBegin;
22 
23   ierr=VecGetOwnershipRange(w,&low,&high);CHKERRQ(ierr);
24   ierr=VecGetOwnershipRange(x,&low1,&high1);CHKERRQ(ierr);
25   ierr=VecGetOwnershipRange(gv,&low2,&high2);CHKERRQ(ierr);
26 
27   *gdx=0.0;
28   ierr = VecGetArray(x,&xptr);CHKERRQ(ierr);
29   ierr = VecGetArray(w,&wptr);CHKERRQ(ierr);
30   ierr = VecGetArray(gv,&gptr);CHKERRQ(ierr);
31 
32   for (i=0;i<high-low;i++)
33     {
34       if (xptr[i]*wptr[i]<0.0)
35 	wptr[i]=0.0;
36 
37       *gdx = *gdx + gptr[i]*(wptr[i]-xptr[i]);
38     }
39   ierr = VecRestoreArray(w,&wptr);CHKERRQ(ierr);
40   ierr = VecRestoreArray(x,&xptr);CHKERRQ(ierr);
41   ierr = VecRestoreArray(gv,&gptr);CHKERRQ(ierr);
42 
43   PetscFunctionReturn(0);
44 }
45 
46 #undef __FUNCT__
47 #define __FUNCT__ "TaoLineSearchDestroy_OWArmijo"
48 static PetscErrorCode TaoLineSearchDestroy_OWArmijo(TaoLineSearch ls)
49 {
50   TAOLINESEARCH_OWARMIJO_CTX *armP = (TAOLINESEARCH_OWARMIJO_CTX *)ls->data;
51   PetscErrorCode ierr;
52 
53   PetscFunctionBegin;
54   ierr = PetscFree(armP->memory); CHKERRQ(ierr);
55   if (armP->x) {
56     ierr = PetscObjectDereference((PetscObject)armP->x); CHKERRQ(ierr);
57   }
58   ierr = VecDestroy(&armP->work); CHKERRQ(ierr);
59   ierr = PetscFree(ls->data); CHKERRQ(ierr);
60   PetscFunctionReturn(0);
61 }
62 
63 #undef __FUNCT__
64 #define __FUNCT__ "TaoLineSearchSetFromOptions_OWArmijo"
65 static PetscErrorCode TaoLineSearchSetFromOptions_OWArmijo(TaoLineSearch ls)
66 {
67   TAOLINESEARCH_OWARMIJO_CTX *armP = (TAOLINESEARCH_OWARMIJO_CTX *)ls->data;
68   PetscErrorCode ierr;
69 
70   PetscFunctionBegin;
71   ierr = PetscOptionsHead("OWArmijo linesearch options");CHKERRQ(ierr);
72   ierr = PetscOptionsReal("-tao_ls_OWArmijo_alpha", "initial reference constant", "", armP->alpha, &armP->alpha, 0); CHKERRQ(ierr);
73   ierr = PetscOptionsReal("-tao_ls_OWArmijo_beta_inf", "decrease constant one", "", armP->beta_inf, &armP->beta_inf, 0); CHKERRQ(ierr);
74   ierr = PetscOptionsReal("-tao_ls_OWArmijo_beta", "decrease constant", "", armP->beta, &armP->beta, 0); CHKERRQ(ierr);
75   ierr = PetscOptionsReal("-tao_ls_OWArmijo_sigma", "acceptance constant", "", armP->sigma, &armP->sigma, 0); CHKERRQ(ierr);
76   ierr = PetscOptionsInt("-tao_ls_OWArmijo_memory_size", "number of historical elements", "", armP->memorySize, &armP->memorySize, 0); CHKERRQ(ierr);
77   ierr = PetscOptionsInt("-tao_ls_OWArmijo_reference_policy", "policy for updating reference value", "", armP->referencePolicy, &armP->referencePolicy, 0); CHKERRQ(ierr);
78   ierr = PetscOptionsInt("-tao_ls_OWArmijo_replacement_policy", "policy for updating memory", "", armP->replacementPolicy, &armP->replacementPolicy, 0); CHKERRQ(ierr);
79   ierr = PetscOptionsBool("-tao_ls_OWArmijo_nondescending","Use nondescending OWArmijo algorithm","",armP->nondescending,&armP->nondescending, 0); CHKERRQ(ierr);
80   ierr = PetscOptionsTail();CHKERRQ(ierr);
81   PetscFunctionReturn(0);
82 }
83 
84 #undef __FUNCT__
85 #define __FUNCT__ "TaoLineSearchView_OWArmijo"
86 static PetscErrorCode TaoLineSearchView_OWArmijo(TaoLineSearch ls, PetscViewer pv)
87 {
88   TAOLINESEARCH_OWARMIJO_CTX *armP = (TAOLINESEARCH_OWARMIJO_CTX *)ls->data;
89   PetscBool isascii;
90   PetscErrorCode ierr;
91 
92   PetscFunctionBegin;
93   ierr = PetscObjectTypeCompare((PetscObject)pv, PETSCVIEWERASCII, &isascii); CHKERRQ(ierr);
94   if (isascii) {
95     ierr = PetscViewerASCIIPrintf(pv,"  maxf=%d, ftol=%g, gtol=%g\n",ls->max_funcs, ls->rtol, ls->ftol); CHKERRQ(ierr);
96     ierr=PetscViewerASCIIPrintf(pv,"  OWArmijo linesearch",armP->alpha);CHKERRQ(ierr);
97     if (armP->nondescending) {
98       ierr = PetscViewerASCIIPrintf(pv, " (nondescending)"); CHKERRQ(ierr);
99     }
100     ierr=PetscViewerASCIIPrintf(pv,": alpha=%g beta=%g ",armP->alpha,armP->beta);CHKERRQ(ierr);
101     ierr=PetscViewerASCIIPrintf(pv,"sigma=%g ",armP->sigma);CHKERRQ(ierr);
102     ierr=PetscViewerASCIIPrintf(pv,"memsize=%d\n",armP->memorySize);CHKERRQ(ierr);
103   } else {
104     SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Viewer type %s not supported for OWArmijo TaoLineSearch",((PetscObject)pv)->type_name);
105   }
106   PetscFunctionReturn(0);
107 }
108 
109 #undef __FUNCT__
110 #define __FUNCT__ "TaoLineSearchApply_OWArmijo"
111 /* @ TaoApply_OWArmijo - This routine performs a linesearch. It
112    backtracks until the (nonmonotone) OWArmijo conditions are satisfied.
113 
114    Input Parameters:
115 +  tao - TAO_SOLVER context
116 .  X - current iterate (on output X contains new iterate, X + step*S)
117 .  S - search direction
118 .  f - merit function evaluated at X
119 .  G - gradient of merit function evaluated at X
120 .  W - work vector
121 -  step - initial estimate of step length
122 
123    Output parameters:
124 +  f - merit function evaluated at new iterate, X + step*S
125 .  G - gradient of merit function evaluated at new iterate, X + step*S
126 .  X - new iterate
127 -  step - final step length
128 
129    Info is set to one of:
130 .   0 - the line search succeeds; the sufficient decrease
131    condition and the directional derivative condition hold
132 
133    negative number if an input parameter is invalid
134 -   -1 -  step < 0
135 
136    positive number > 1 if the line search otherwise terminates
137 +    1 -  Step is at the lower bound, stepmin.
138 @ */
139 static PetscErrorCode TaoLineSearchApply_OWArmijo(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, Vec s)
140 {
141   TAOLINESEARCH_OWARMIJO_CTX *armP = (TAOLINESEARCH_OWARMIJO_CTX *)ls->data;
142   PetscErrorCode ierr;
143   PetscInt i;
144   PetscReal fact, ref, gdx;
145   PetscInt idx;
146   PetscBool g_computed=PETSC_FALSE; /* to prevent extra gradient computation */
147 
148   Vec g_old;
149   PetscReal owlqn_minstep=0.005;
150   PetscReal partgdx;
151 
152   PetscFunctionBegin;
153 
154   fact = 0.0;
155   ls->nfeval=0;
156   ls->reason = TAOLINESEARCH_CONTINUE_ITERATING;
157   if (!armP->work) {
158     ierr = VecDuplicate(x,&armP->work); CHKERRQ(ierr);
159     armP->x = x;
160     ierr = PetscObjectReference((PetscObject)armP->x); CHKERRQ(ierr);
161   }
162   /* If x has changed, then recreate work */
163   else if (x != armP->x) {
164     ierr = VecDestroy(&armP->work); CHKERRQ(ierr);
165     ierr = VecDuplicate(x,&armP->work); CHKERRQ(ierr);
166     ierr = PetscObjectDereference((PetscObject)armP->x); CHKERRQ(ierr);
167     armP->x = x;
168     ierr = PetscObjectReference((PetscObject)armP->x); CHKERRQ(ierr);
169   }
170 
171   /* Check linesearch parameters */
172   if (armP->alpha < 1) {
173     ierr = PetscInfo1(ls,"OWArmijo line search error: alpha (%g) < 1\n", armP->alpha); CHKERRQ(ierr);
174     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
175   }
176 
177   else if ((armP->beta <= 0) || (armP->beta >= 1)) {
178     ierr = PetscInfo1(ls,"OWArmijo line search error: beta (%g) invalid\n", armP->beta); CHKERRQ(ierr);
179     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
180 
181   }
182 
183   else if ((armP->beta_inf <= 0) || (armP->beta_inf >= 1)) {
184     ierr = PetscInfo1(ls,"OWArmijo line search error: beta_inf (%g) invalid\n", armP->beta_inf); CHKERRQ(ierr);
185     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
186   }
187 
188   else if ((armP->sigma <= 0) || (armP->sigma >= 0.5)) {
189     ierr = PetscInfo1(ls,"OWArmijo line search error: sigma (%g) invalid\n", armP->sigma); CHKERRQ(ierr);
190     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
191   }
192 
193   else if (armP->memorySize < 1) {
194     ierr = PetscInfo1(ls,"OWArmijo line search error: memory_size (%d) < 1\n", armP->memorySize); CHKERRQ(ierr);
195     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
196   }
197 
198   else if ((armP->referencePolicy != REFERENCE_MAX) &&
199       (armP->referencePolicy != REFERENCE_AVE) &&
200       (armP->referencePolicy != REFERENCE_MEAN)) {
201     ierr = PetscInfo(ls,"OWArmijo line search error: reference_policy invalid\n"); CHKERRQ(ierr);
202     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
203 
204   }
205 
206   else if ((armP->replacementPolicy != REPLACE_FIFO) &&
207       (armP->replacementPolicy != REPLACE_MRU)) {
208     ierr = PetscInfo(ls,"OWArmijo line search error: replacement_policy invalid\n"); CHKERRQ(ierr);
209     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
210   }
211 
212   else if (PetscIsInfOrNanReal(*f)) {
213     ierr = PetscInfo(ls,"OWArmijo line search error: initial function inf or nan\n"); CHKERRQ(ierr);
214     ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
215   }
216 
217   if (ls->reason != TAOLINESEARCH_CONTINUE_ITERATING) {
218     PetscFunctionReturn(0);
219   }
220 
221   /* Check to see of the memory has been allocated.  If not, allocate
222      the historical array and populate it with the initial function
223      values. */
224   if (!armP->memory) {
225     ierr = PetscMalloc(sizeof(PetscReal)*armP->memorySize, &armP->memory ); CHKERRQ(ierr);
226   }
227 
228   if (!armP->memorySetup) {
229     for (i = 0; i < armP->memorySize; i++) {
230       armP->memory[i] = armP->alpha*(*f);
231     }
232 
233     armP->current = 0;
234     armP->lastReference = armP->memory[0];
235     armP->memorySetup=PETSC_TRUE;
236   }
237 
238   /* Calculate reference value (MAX) */
239   ref = armP->memory[0];
240   idx = 0;
241 
242   for (i = 1; i < armP->memorySize; i++) {
243     if (armP->memory[i] > ref) {
244       ref = armP->memory[i];
245       idx = i;
246     }
247   }
248 
249   if (armP->referencePolicy == REFERENCE_AVE) {
250     ref = 0;
251     for (i = 0; i < armP->memorySize; i++) {
252       ref += armP->memory[i];
253     }
254     ref = ref / armP->memorySize;
255     ref = PetscMax(ref, armP->memory[armP->current]);
256   }
257   else if (armP->referencePolicy == REFERENCE_MEAN) {
258     ref = PetscMin(ref, 0.5*(armP->lastReference + armP->memory[armP->current]));
259   }
260 
261 
262   if (armP->nondescending) {
263     fact = armP->sigma;
264   }
265 
266   VecDuplicate(g,&g_old);
267   VecCopy(g,g_old);
268 
269   ls->step = ls->initstep;
270 
271   while (ls->step >= owlqn_minstep && ls->nfeval < ls->max_funcs) {
272     /* Calculate iterate */
273     ierr = VecCopy(x,armP->work); CHKERRQ(ierr);
274     ierr = VecAXPY(armP->work,ls->step,s); CHKERRQ(ierr);
275 
276     partgdx=0.0;
277     ierr = ProjWork_OWLQN(armP->work,x,g_old,&partgdx);
278     ierr = (PetscErrorCode)MPI_Allreduce(&partgdx,&gdx,1,MPI_DOUBLE,MPI_SUM,PETSC_COMM_WORLD); CHKERRQ(ierr);
279 
280     /* Check the condition of gdx */
281     if (PetscIsInfOrNanReal(gdx)) {
282       ierr = PetscInfo1(ls,"Initial Line Search step * g is Inf or Nan (%g)\n",gdx); CHKERRQ(ierr);
283       ls->reason=TAOLINESEARCH_FAILED_INFORNAN;
284       PetscFunctionReturn(0);
285     }
286     if (gdx >= 0.0) {
287       ierr = PetscInfo1(ls,"Initial Line Search step is not descent direction (g's=%g)\n",gdx); CHKERRQ(ierr);
288       ls->reason = TAOLINESEARCH_FAILED_ASCENT;
289       PetscFunctionReturn(0);
290     }
291 
292     /* Calculate function at new iterate */
293     ierr = TaoLineSearchComputeObjectiveAndGradient(ls,armP->work,f,g); CHKERRQ(ierr);
294     g_computed=PETSC_TRUE;
295 
296 
297     if (ls->step == ls->initstep) {
298       ls->f_fullstep = *f;
299     }
300 
301     if (PetscIsInfOrNanReal(*f)) {
302       ls->step *= armP->beta_inf;
303     }
304     else {
305       /* Check descent condition */
306       if (armP->nondescending && *f <= ref - ls->step*fact*ref)
307 	break;
308       if (!armP->nondescending && *f <= ref + armP->sigma * gdx) {
309         break;
310       }
311 
312       ls->step *= armP->beta;
313     }
314   }
315 
316   VecDestroy(&g_old);
317 
318   /* Check termination */
319   if (PetscIsInfOrNanReal(*f)) {
320     ierr = PetscInfo(ls, "Function is inf or nan.\n"); CHKERRQ(ierr);
321     ls->reason = TAOLINESEARCH_FAILED_BADPARAMETER;
322   } else if (ls->step < owlqn_minstep) {
323     ierr = PetscInfo(ls, "Step length is below tolerance.\n"); CHKERRQ(ierr);
324     ls->reason = TAOLINESEARCH_HALTED_RTOL;
325   } else if (ls->nfeval >= ls->max_funcs) {
326     ierr = PetscInfo2(ls, "Number of line search function evals (%d) > maximum allowed (%d)\n",ls->nfeval, ls->max_funcs); CHKERRQ(ierr);
327     ls->reason = TAOLINESEARCH_HALTED_MAXFCN;
328   }
329   if (ls->reason) {
330     PetscFunctionReturn(0);
331   }
332 
333   /* Successful termination, update memory */
334   armP->lastReference = ref;
335   if (armP->replacementPolicy == REPLACE_FIFO) {
336     armP->memory[armP->current++] = *f;
337     if (armP->current >= armP->memorySize) {
338       armP->current = 0;
339     }
340   } else {
341     armP->current = idx;
342     armP->memory[idx] = *f;
343   }
344 
345   /* Update iterate and compute gradient */
346   ierr = VecCopy(armP->work,x); CHKERRQ(ierr);
347   if (!g_computed) {
348     ierr = TaoLineSearchComputeGradient(ls, x, g); CHKERRQ(ierr);
349   }
350 
351   /* Finish computations */
352   ierr = PetscInfo2(ls, "%d function evals in line search, step = %10.4f\n",ls->nfeval, ls->step); CHKERRQ(ierr);
353 
354   PetscFunctionReturn(0);
355 }
356 
357 EXTERN_C_BEGIN
358 #undef __FUNCT__
359 #define __FUNCT__ "TaoLineSearchCreate_OWArmijo"
360 PetscErrorCode TaoLineSearchCreate_OWArmijo(TaoLineSearch ls)
361 {
362   TAOLINESEARCH_OWARMIJO_CTX *armP;
363   PetscErrorCode ierr;
364 
365   PetscFunctionBegin;
366   PetscValidHeaderSpecific(ls,TAOLINESEARCH_CLASSID,1);
367   ierr = PetscNewLog(ls,&armP);CHKERRQ(ierr);
368 
369   armP->memory = NULL;
370   armP->alpha = 1.0;
371   armP->beta = 0.25;
372   armP->beta_inf = 0.25;
373   armP->sigma = 1e-4;
374   armP->memorySize = 1;
375   armP->referencePolicy = REFERENCE_MAX;
376   armP->replacementPolicy = REPLACE_MRU;
377   armP->nondescending=PETSC_FALSE;
378   ls->data = (void*)armP;
379   ls->initstep=0.1;
380   ls->ops->setup=0;
381   ls->ops->apply=TaoLineSearchApply_OWArmijo;
382   ls->ops->view = TaoLineSearchView_OWArmijo;
383   ls->ops->destroy = TaoLineSearchDestroy_OWArmijo;
384   ls->ops->setfromoptions = TaoLineSearchSetFromOptions_OWArmijo;
385 
386   PetscFunctionReturn(0);
387 }
388 EXTERN_C_END
389