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: }