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