Actual source code: armijo.c
1: #include "petscvec.h"
2: #include "taosolver.h"
3: #include "private/taolinesearch_impl.h"
4: #include "armijo.h"
6: #define REPLACE_FIFO 1
7: #define REPLACE_MRU 2
9: #define REFERENCE_MAX 1
10: #define REFERENCE_AVE 2
11: #define REFERENCE_MEAN 3
15: static PetscErrorCode TaoLineSearchDestroy_Armijo(TaoLineSearch ls)
16: {
17: TAOLINESEARCH_ARMIJO_CTX *armP = (TAOLINESEARCH_ARMIJO_CTX *)ls->data;
22: if (armP->memory != PETSC_NULL) {
23: PetscFree(armP->memory);
24: armP->memory = PETSC_NULL;
25: }
26: if (armP->x) {
27: PetscObjectDereference((PetscObject)armP->x);
28: }
29: VecDestroy(&armP->work);
30: PetscFree(ls->data);
31: ls->data = PETSC_NULL;
32: return(0);
33: }
37: static PetscErrorCode TaoLineSearchReset_Armijo(TaoLineSearch ls)
38: {
39: TAOLINESEARCH_ARMIJO_CTX *armP = (TAOLINESEARCH_ARMIJO_CTX *)ls->data;
44: if (armP->memory != PETSC_NULL) {
45: PetscFree(armP->memory);
46: armP->memory = PETSC_NULL;
47: }
48: armP->memorySetup = PETSC_FALSE;
49: return(0);
50: }
54: static PetscErrorCode TaoLineSearchSetFromOptions_Armijo(TaoLineSearch ls)
55: {
56: TAOLINESEARCH_ARMIJO_CTX *armP = (TAOLINESEARCH_ARMIJO_CTX *)ls->data;
60: PetscOptionsHead("Armijo linesearch options");
61: PetscOptionsReal("-tao_ls_armijo_alpha", "initial reference constant", "", armP->alpha, &armP->alpha, 0);
62: PetscOptionsReal("-tao_ls_armijo_beta_inf", "decrease constant one", "", armP->beta_inf, &armP->beta_inf, 0);
63: PetscOptionsReal("-tao_ls_armijo_beta", "decrease constant", "", armP->beta, &armP->beta, 0);
64: PetscOptionsReal("-tao_ls_armijo_sigma", "acceptance constant", "", armP->sigma, &armP->sigma, 0);
65: PetscOptionsInt("-tao_ls_armijo_memory_size", "number of historical elements", "", armP->memorySize, &armP->memorySize, 0);
66: PetscOptionsInt("-tao_ls_armijo_reference_policy", "policy for updating reference value", "", armP->referencePolicy, &armP->referencePolicy, 0);
67: PetscOptionsInt("-tao_ls_armijo_replacement_policy", "policy for updating memory", "", armP->replacementPolicy, &armP->replacementPolicy, 0);
68: PetscOptionsBool("-tao_ls_armijo_nondescending","Use nondescending armijo algorithm","",armP->nondescending,&armP->nondescending, 0);
69: PetscOptionsTail();
70: return(0);
71: }
75: static PetscErrorCode TaoLineSearchView_Armijo(TaoLineSearch ls, PetscViewer pv)
76: {
77: TAOLINESEARCH_ARMIJO_CTX *armP = (TAOLINESEARCH_ARMIJO_CTX *)ls->data;
78: PetscBool isascii;
82: PetscTypeCompare((PetscObject)pv, PETSCVIEWERASCII, &isascii);
83: if (isascii) {
84: PetscViewerASCIIPrintf(pv," maxf=%D, ftol=%G, gtol=%G\n",ls->max_funcs, ls->rtol, ls->ftol);
85: ierr=PetscViewerASCIIPrintf(pv," Armijo linesearch",armP->alpha);
86: if (armP->nondescending) {
87: PetscViewerASCIIPrintf(pv, " (nondescending)");
88: }
89: if (ls->bounded) {
90: PetscViewerASCIIPrintf(pv," (projected)");
91: }
92: ierr=PetscViewerASCIIPrintf(pv,": alpha=%G beta=%G ",armP->alpha,armP->beta);
93: ierr=PetscViewerASCIIPrintf(pv,"sigma=%G ",armP->sigma);
94: ierr=PetscViewerASCIIPrintf(pv,"memsize=%D\n",armP->memorySize);
95: } else {
96: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Viewer type %s not supported for Armijo TaoLineSearch",((PetscObject)pv)->type_name);
97: }
98: return(0);
99: }
103: /* @ TaoApply_Armijo - This routine performs a linesearch. It
104: backtracks until the (nonmonotone) Armijo conditions are satisfied.
106: Input Parameters:
107: + tao - TaoSolver context
108: . X - current iterate (on output X contains new iterate, X + step*S)
109: . S - search direction
110: . f - merit function evaluated at X
111: . G - gradient of merit function evaluated at X
112: . W - work vector
113: - step - initial estimate of step length
115: Output parameters:
116: + f - merit function evaluated at new iterate, X + step*S
117: . G - gradient of merit function evaluated at new iterate, X + step*S
118: . X - new iterate
119: - step - final step length
121: Info is set to one of:
122: . 0 - the line search succeeds; the sufficient decrease
123: condition and the directional derivative condition hold
125: negative number if an input parameter is invalid
126: - -1 - step < 0
128: positive number > 1 if the line search otherwise terminates
129: + 1 - Step is at the lower bound, stepmin.
130: @ */
131: static PetscErrorCode TaoLineSearchApply_Armijo(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, Vec s)
132: {
133: TAOLINESEARCH_ARMIJO_CTX *armP = (TAOLINESEARCH_ARMIJO_CTX *)ls->data;
135: PetscInt i;
136: PetscReal fact, ref, gdx;
137: PetscInt idx;
138: PetscBool g_computed=PETSC_FALSE; /* to prevent extra gradient computation */
142: ls->reason = TAOLINESEARCH_CONTINUE_ITERATING;
143: if (!armP->work) {
144: VecDuplicate(x,&armP->work);
145: armP->x = x;
146: PetscObjectReference((PetscObject)armP->x);
147: }
148: /* If x has changed, then recreate work */
149: else if (x != armP->x) {
150: VecDestroy(&armP->work);
151: VecDuplicate(x,&armP->work);
152: PetscObjectDereference((PetscObject)armP->x);
153: armP->x = x;
154: PetscObjectReference((PetscObject)armP->x);
155: }
157: /* Check linesearch parameters */
158: if (armP->alpha < 1) {
159: PetscInfo1(ls,"Armijo line search error: alpha (%G) < 1\n", armP->alpha);
160: ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
161: }
162:
163: else if ((armP->beta <= 0) || (armP->beta >= 1)) {
164: PetscInfo1(ls,"Armijo line search error: beta (%G) invalid\n", armP->beta);
165: ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
167: }
168:
169: else if ((armP->beta_inf <= 0) || (armP->beta_inf >= 1)) {
170: PetscInfo1(ls,"Armijo line search error: beta_inf (%G) invalid\n", armP->beta_inf);
171: ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
172: }
174: else if ((armP->sigma <= 0) || (armP->sigma >= 0.5)) {
175: PetscInfo1(ls,"Armijo line search error: sigma (%G) invalid\n", armP->sigma);
176: ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
177: }
178:
179: else if (armP->memorySize < 1) {
180: PetscInfo1(ls,"Armijo line search error: memory_size (%D) < 1\n", armP->memorySize);
181: ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
182: }
183:
184: else if ((armP->referencePolicy != REFERENCE_MAX) &&
185: (armP->referencePolicy != REFERENCE_AVE) &&
186: (armP->referencePolicy != REFERENCE_MEAN)) {
187: PetscInfo(ls,"Armijo line search error: reference_policy invalid\n");
188: ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
190: }
191:
192: else if ((armP->replacementPolicy != REPLACE_FIFO) &&
193: (armP->replacementPolicy != REPLACE_MRU)) {
194: PetscInfo(ls,"Armijo line search error: replacement_policy invalid\n");
195: ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
196: }
197:
198: else if (PetscIsInfOrNanReal(*f)) {
199: PetscInfo(ls,"Armijo line search error: initial function inf or nan\n");
200: ls->reason=TAOLINESEARCH_FAILED_BADPARAMETER;
201: }
203: if (ls->reason != TAOLINESEARCH_CONTINUE_ITERATING) {
204: return(0);
205: }
207: /* Check to see of the memory has been allocated. If not, allocate
208: the historical array and populate it with the initial function
209: values. */
210: if (armP->memory == PETSC_NULL) {
211: PetscMalloc(sizeof(PetscReal)*armP->memorySize, &armP->memory );
212: }
214: if (!armP->memorySetup) {
215: for (i = 0; i < armP->memorySize; i++) {
216: armP->memory[i] = armP->alpha*(*f);
217: }
219: armP->current = 0;
220: armP->lastReference = armP->memory[0];
221: armP->memorySetup=PETSC_TRUE;
222: }
223:
224: /* Calculate reference value (MAX) */
225: ref = armP->memory[0];
226: idx = 0;
228: for (i = 1; i < armP->memorySize; i++) {
229: if (armP->memory[i] > ref) {
230: ref = armP->memory[i];
231: idx = i;
232: }
233: }
235: if (armP->referencePolicy == REFERENCE_AVE) {
236: ref = 0;
237: for (i = 0; i < armP->memorySize; i++) {
238: ref += armP->memory[i];
239: }
240: ref = ref / armP->memorySize;
241: ref = PetscMax(ref, armP->memory[armP->current]);
242: }
243: else if (armP->referencePolicy == REFERENCE_MEAN) {
244: ref = PetscMin(ref, 0.5*(armP->lastReference + armP->memory[armP->current]));
245: }
246: VecDot(g,s,&gdx);
248: if (PetscIsInfOrNanReal(gdx)) {
249: PetscInfo1(ls,"Initial Line Search step * g is Inf or Nan (%G)\n",gdx);
250: ls->reason=TAOLINESEARCH_FAILED_INFORNAN;
251: return(0);
252: }
253: if (gdx >= 0.0) {
254: PetscInfo1(ls,"Initial Line Search step is not descent direction (g's=%G)\n",gdx);
255: ls->reason = TAOLINESEARCH_FAILED_ASCENT;
256: return(0);
257: }
258:
259: if (armP->nondescending) {
260: fact = armP->sigma;
261: } else {
262: fact = armP->sigma * gdx;
263: }
264: ls->step = ls->initstep;
265: while (ls->step >= ls->stepmin && (ls->nfeval+ls->nfgeval) < ls->max_funcs) {
266: /* Calculate iterate */
267: VecCopy(x,armP->work);
268: VecAXPY(armP->work,ls->step,s);
269: if (ls->bounded) {
270: VecMedian(ls->lower,armP->work,ls->upper,armP->work);
271: }
273: /* Calculate function at new iterate */
274: if (ls->hasobjective) {
275: TaoLineSearchComputeObjective(ls,armP->work,f);
276: g_computed=PETSC_FALSE;
277: } else if (ls->usegts) {
278: TaoLineSearchComputeObjectiveAndGTS(ls,armP->work,f,&gdx);
279: g_computed=PETSC_FALSE;
280: } else {
281: TaoLineSearchComputeObjectiveAndGradient(ls,armP->work,f,g);
282: g_computed=PETSC_TRUE;
283: }
284: if (ls->step == ls->initstep) {
285: ls->f_fullstep = *f;
286: }
288: if (PetscIsInfOrNanReal(*f)) {
289: ls->step *= armP->beta_inf;
290: }
291: else {
292: /* Check descent condition */
293: if (armP->nondescending && *f <= ref - ls->step*fact*ref)
294: break;
295: if (!armP->nondescending && *f <= ref + ls->step*fact) {
296: break;
297: }
299: ls->step *= armP->beta;
300: }
301: }
303: /* Check termination */
304: if (PetscIsInfOrNanReal(*f)) {
305: PetscInfo(ls, "Function is inf or nan.\n");
306: ls->reason = TAOLINESEARCH_FAILED_INFORNAN;
307: } else if (ls->step < ls->stepmin) {
308: PetscInfo(ls, "Step length is below tolerance.\n");
309: ls->reason = TAOLINESEARCH_HALTED_RTOL;
310: } else if ((ls->nfeval+ls->nfgeval) >= ls->max_funcs) {
311: PetscInfo2(ls, "Number of line search function evals (%D) > maximum allowed (%D)\n",ls->nfeval+ls->nfgeval, ls->max_funcs);
312: ls->reason = TAOLINESEARCH_HALTED_MAXFCN;
313: }
314: if (ls->reason) {
315: return(0);
316: }
318: /* Successful termination, update memory */
319: armP->lastReference = ref;
320: if (armP->replacementPolicy == REPLACE_FIFO) {
321: armP->memory[armP->current++] = *f;
322: if (armP->current >= armP->memorySize) {
323: armP->current = 0;
324: }
325: } else {
326: armP->current = idx;
327: armP->memory[idx] = *f;
328: }
330: /* Update iterate and compute gradient */
331: VecCopy(armP->work,x);
332: if (!g_computed) {
333: TaoLineSearchComputeGradient(ls, x, g);
334: }
336: /* Finish computations */
337: PetscInfo2(ls, "%D function evals in line search, step = %G\n",ls->nfeval, ls->step);
339: return(0);
340: }
345: PetscErrorCode TaoLineSearchCreate_Armijo(TaoLineSearch ls)
346: {
347: TAOLINESEARCH_ARMIJO_CTX *armP;
352: PetscNewLog(ls,TAOLINESEARCH_ARMIJO_CTX, &armP);
354: armP->memory = PETSC_NULL;
355: armP->alpha = 1.0;
356: armP->beta = 0.5;
357: armP->beta_inf = 0.5;
358: armP->sigma = 1e-4;
359: armP->memorySize = 1;
360: armP->referencePolicy = REFERENCE_MAX;
361: armP->replacementPolicy = REPLACE_MRU;
362: armP->nondescending=PETSC_FALSE;
363: ls->data = (void*)armP;
364: ls->initstep=1.0;
365: ls->ops->setup=0;
366: ls->ops->apply=TaoLineSearchApply_Armijo;
367: ls->ops->view = TaoLineSearchView_Armijo;
368: ls->ops->destroy = TaoLineSearchDestroy_Armijo;
369: ls->ops->reset = TaoLineSearchReset_Armijo;
370: ls->ops->setfromoptions = TaoLineSearchSetFromOptions_Armijo;
372: return(0);
373: }