1 2 #include <petsc/private/taolinesearchimpl.h> 3 #include <../src/tao/linesearch/impls/owarmijo/owarmijo.h> 4 5 #define REPLACE_FIFO 1 6 #define REPLACE_MRU 2 7 8 #define REFERENCE_MAX 1 9 #define REFERENCE_AVE 2 10 #define REFERENCE_MEAN 3 11 12 static PetscErrorCode ProjWork_OWLQN(Vec w, Vec x, Vec gv, PetscReal *gdx) { 13 const PetscReal *xptr, *gptr; 14 PetscReal *wptr; 15 PetscInt low, high, low1, high1, low2, high2, i; 16 17 PetscFunctionBegin; 18 PetscCall(VecGetOwnershipRange(w, &low, &high)); 19 PetscCall(VecGetOwnershipRange(x, &low1, &high1)); 20 PetscCall(VecGetOwnershipRange(gv, &low2, &high2)); 21 22 *gdx = 0.0; 23 PetscCall(VecGetArray(w, &wptr)); 24 PetscCall(VecGetArrayRead(x, &xptr)); 25 PetscCall(VecGetArrayRead(gv, &gptr)); 26 27 for (i = 0; i < high - low; i++) { 28 if (xptr[i] * wptr[i] < 0.0) wptr[i] = 0.0; 29 *gdx = *gdx + gptr[i] * (wptr[i] - xptr[i]); 30 } 31 PetscCall(VecRestoreArray(w, &wptr)); 32 PetscCall(VecRestoreArrayRead(x, &xptr)); 33 PetscCall(VecRestoreArrayRead(gv, &gptr)); 34 PetscFunctionReturn(0); 35 } 36 37 static PetscErrorCode TaoLineSearchDestroy_OWArmijo(TaoLineSearch ls) { 38 TaoLineSearch_OWARMIJO *armP = (TaoLineSearch_OWARMIJO *)ls->data; 39 40 PetscFunctionBegin; 41 PetscCall(PetscFree(armP->memory)); 42 if (armP->x) PetscCall(PetscObjectDereference((PetscObject)armP->x)); 43 PetscCall(VecDestroy(&armP->work)); 44 PetscCall(PetscFree(ls->data)); 45 PetscFunctionReturn(0); 46 } 47 48 static PetscErrorCode TaoLineSearchSetFromOptions_OWArmijo(TaoLineSearch ls, PetscOptionItems *PetscOptionsObject) { 49 TaoLineSearch_OWARMIJO *armP = (TaoLineSearch_OWARMIJO *)ls->data; 50 51 PetscFunctionBegin; 52 PetscOptionsHeadBegin(PetscOptionsObject, "OWArmijo linesearch options"); 53 PetscCall(PetscOptionsReal("-tao_ls_OWArmijo_alpha", "initial reference constant", "", armP->alpha, &armP->alpha, NULL)); 54 PetscCall(PetscOptionsReal("-tao_ls_OWArmijo_beta_inf", "decrease constant one", "", armP->beta_inf, &armP->beta_inf, NULL)); 55 PetscCall(PetscOptionsReal("-tao_ls_OWArmijo_beta", "decrease constant", "", armP->beta, &armP->beta, NULL)); 56 PetscCall(PetscOptionsReal("-tao_ls_OWArmijo_sigma", "acceptance constant", "", armP->sigma, &armP->sigma, NULL)); 57 PetscCall(PetscOptionsInt("-tao_ls_OWArmijo_memory_size", "number of historical elements", "", armP->memorySize, &armP->memorySize, NULL)); 58 PetscCall(PetscOptionsInt("-tao_ls_OWArmijo_reference_policy", "policy for updating reference value", "", armP->referencePolicy, &armP->referencePolicy, NULL)); 59 PetscCall(PetscOptionsInt("-tao_ls_OWArmijo_replacement_policy", "policy for updating memory", "", armP->replacementPolicy, &armP->replacementPolicy, NULL)); 60 PetscCall(PetscOptionsBool("-tao_ls_OWArmijo_nondescending", "Use nondescending OWArmijo algorithm", "", armP->nondescending, &armP->nondescending, NULL)); 61 PetscOptionsHeadEnd(); 62 PetscFunctionReturn(0); 63 } 64 65 static PetscErrorCode TaoLineSearchView_OWArmijo(TaoLineSearch ls, PetscViewer pv) { 66 TaoLineSearch_OWARMIJO *armP = (TaoLineSearch_OWARMIJO *)ls->data; 67 PetscBool isascii; 68 69 PetscFunctionBegin; 70 PetscCall(PetscObjectTypeCompare((PetscObject)pv, PETSCVIEWERASCII, &isascii)); 71 if (isascii) { 72 PetscCall(PetscViewerASCIIPrintf(pv, " OWArmijo linesearch")); 73 if (armP->nondescending) { PetscCall(PetscViewerASCIIPrintf(pv, " (nondescending)")); } 74 PetscCall(PetscViewerASCIIPrintf(pv, ": alpha=%g beta=%g ", (double)armP->alpha, (double)armP->beta)); 75 PetscCall(PetscViewerASCIIPrintf(pv, "sigma=%g ", (double)armP->sigma)); 76 PetscCall(PetscViewerASCIIPrintf(pv, "memsize=%" PetscInt_FMT "\n", armP->memorySize)); 77 } 78 PetscFunctionReturn(0); 79 } 80 81 /* @ TaoApply_OWArmijo - This routine performs a linesearch. It 82 backtracks until the (nonmonotone) OWArmijo conditions are satisfied. 83 84 Input Parameters: 85 + tao - TAO_SOLVER context 86 . X - current iterate (on output X contains new iterate, X + step*S) 87 . S - search direction 88 . f - merit function evaluated at X 89 . G - gradient of merit function evaluated at X 90 . W - work vector 91 - step - initial estimate of step length 92 93 Output parameters: 94 + f - merit function evaluated at new iterate, X + step*S 95 . G - gradient of merit function evaluated at new iterate, X + step*S 96 . X - new iterate 97 - step - final step length 98 99 Info is set to one of: 100 . 0 - the line search succeeds; the sufficient decrease 101 condition and the directional derivative condition hold 102 103 negative number if an input parameter is invalid 104 - -1 - step < 0 105 106 positive number > 1 if the line search otherwise terminates 107 + 1 - Step is at the lower bound, stepmin. 108 @ */ 109 static PetscErrorCode TaoLineSearchApply_OWArmijo(TaoLineSearch ls, Vec x, PetscReal *f, Vec g, Vec s) { 110 TaoLineSearch_OWARMIJO *armP = (TaoLineSearch_OWARMIJO *)ls->data; 111 PetscInt i, its = 0; 112 PetscReal fact, ref, gdx; 113 PetscInt idx; 114 PetscBool g_computed = PETSC_FALSE; /* to prevent extra gradient computation */ 115 Vec g_old; 116 PetscReal owlqn_minstep = 0.005; 117 PetscReal partgdx; 118 MPI_Comm comm; 119 120 PetscFunctionBegin; 121 PetscCall(PetscObjectGetComm((PetscObject)ls, &comm)); 122 fact = 0.0; 123 ls->nfeval = 0; 124 ls->reason = TAOLINESEARCH_CONTINUE_ITERATING; 125 if (!armP->work) { 126 PetscCall(VecDuplicate(x, &armP->work)); 127 armP->x = x; 128 PetscCall(PetscObjectReference((PetscObject)armP->x)); 129 } else if (x != armP->x) { 130 PetscCall(VecDestroy(&armP->work)); 131 PetscCall(VecDuplicate(x, &armP->work)); 132 PetscCall(PetscObjectDereference((PetscObject)armP->x)); 133 armP->x = x; 134 PetscCall(PetscObjectReference((PetscObject)armP->x)); 135 } 136 137 PetscCall(TaoLineSearchMonitor(ls, 0, *f, 0.0)); 138 139 /* Check linesearch parameters */ 140 if (armP->alpha < 1) { 141 PetscCall(PetscInfo(ls, "OWArmijo line search error: alpha (%g) < 1\n", (double)armP->alpha)); 142 ls->reason = TAOLINESEARCH_FAILED_BADPARAMETER; 143 } else if ((armP->beta <= 0) || (armP->beta >= 1)) { 144 PetscCall(PetscInfo(ls, "OWArmijo line search error: beta (%g) invalid\n", (double)armP->beta)); 145 ls->reason = TAOLINESEARCH_FAILED_BADPARAMETER; 146 } else if ((armP->beta_inf <= 0) || (armP->beta_inf >= 1)) { 147 PetscCall(PetscInfo(ls, "OWArmijo line search error: beta_inf (%g) invalid\n", (double)armP->beta_inf)); 148 ls->reason = TAOLINESEARCH_FAILED_BADPARAMETER; 149 } else if ((armP->sigma <= 0) || (armP->sigma >= 0.5)) { 150 PetscCall(PetscInfo(ls, "OWArmijo line search error: sigma (%g) invalid\n", (double)armP->sigma)); 151 ls->reason = TAOLINESEARCH_FAILED_BADPARAMETER; 152 } else if (armP->memorySize < 1) { 153 PetscCall(PetscInfo(ls, "OWArmijo line search error: memory_size (%" PetscInt_FMT ") < 1\n", armP->memorySize)); 154 ls->reason = TAOLINESEARCH_FAILED_BADPARAMETER; 155 } else if ((armP->referencePolicy != REFERENCE_MAX) && (armP->referencePolicy != REFERENCE_AVE) && (armP->referencePolicy != REFERENCE_MEAN)) { 156 PetscCall(PetscInfo(ls, "OWArmijo line search error: reference_policy invalid\n")); 157 ls->reason = TAOLINESEARCH_FAILED_BADPARAMETER; 158 } else if ((armP->replacementPolicy != REPLACE_FIFO) && (armP->replacementPolicy != REPLACE_MRU)) { 159 PetscCall(PetscInfo(ls, "OWArmijo line search error: replacement_policy invalid\n")); 160 ls->reason = TAOLINESEARCH_FAILED_BADPARAMETER; 161 } else if (PetscIsInfOrNanReal(*f)) { 162 PetscCall(PetscInfo(ls, "OWArmijo line search error: initial function inf or nan\n")); 163 ls->reason = TAOLINESEARCH_FAILED_BADPARAMETER; 164 } 165 166 if (ls->reason != TAOLINESEARCH_CONTINUE_ITERATING) PetscFunctionReturn(0); 167 168 /* Check to see of the memory has been allocated. If not, allocate 169 the historical array and populate it with the initial function 170 values. */ 171 if (!armP->memory) { PetscCall(PetscMalloc1(armP->memorySize, &armP->memory)); } 172 173 if (!armP->memorySetup) { 174 for (i = 0; i < armP->memorySize; i++) { armP->memory[i] = armP->alpha * (*f); } 175 armP->current = 0; 176 armP->lastReference = armP->memory[0]; 177 armP->memorySetup = PETSC_TRUE; 178 } 179 180 /* Calculate reference value (MAX) */ 181 ref = armP->memory[0]; 182 idx = 0; 183 184 for (i = 1; i < armP->memorySize; i++) { 185 if (armP->memory[i] > ref) { 186 ref = armP->memory[i]; 187 idx = i; 188 } 189 } 190 191 if (armP->referencePolicy == REFERENCE_AVE) { 192 ref = 0; 193 for (i = 0; i < armP->memorySize; i++) { ref += armP->memory[i]; } 194 ref = ref / armP->memorySize; 195 ref = PetscMax(ref, armP->memory[armP->current]); 196 } else if (armP->referencePolicy == REFERENCE_MEAN) { 197 ref = PetscMin(ref, 0.5 * (armP->lastReference + armP->memory[armP->current])); 198 } 199 200 if (armP->nondescending) { fact = armP->sigma; } 201 202 PetscCall(VecDuplicate(g, &g_old)); 203 PetscCall(VecCopy(g, g_old)); 204 205 ls->step = ls->initstep; 206 while (ls->step >= owlqn_minstep && ls->nfeval < ls->max_funcs) { 207 /* Calculate iterate */ 208 ++its; 209 PetscCall(VecWAXPY(armP->work, ls->step, s, x)); 210 211 partgdx = 0.0; 212 PetscCall(ProjWork_OWLQN(armP->work, x, g_old, &partgdx)); 213 PetscCall(MPIU_Allreduce(&partgdx, &gdx, 1, MPIU_REAL, MPIU_SUM, comm)); 214 215 /* Check the condition of gdx */ 216 if (PetscIsInfOrNanReal(gdx)) { 217 PetscCall(PetscInfo(ls, "Initial Line Search step * g is Inf or Nan (%g)\n", (double)gdx)); 218 ls->reason = TAOLINESEARCH_FAILED_INFORNAN; 219 PetscFunctionReturn(0); 220 } 221 if (gdx >= 0.0) { 222 PetscCall(PetscInfo(ls, "Initial Line Search step is not descent direction (g's=%g)\n", (double)gdx)); 223 ls->reason = TAOLINESEARCH_FAILED_ASCENT; 224 PetscFunctionReturn(0); 225 } 226 227 /* Calculate function at new iterate */ 228 PetscCall(TaoLineSearchComputeObjectiveAndGradient(ls, armP->work, f, g)); 229 g_computed = PETSC_TRUE; 230 231 PetscCall(TaoLineSearchMonitor(ls, its, *f, ls->step)); 232 233 if (ls->step == ls->initstep) { ls->f_fullstep = *f; } 234 235 if (PetscIsInfOrNanReal(*f)) { 236 ls->step *= armP->beta_inf; 237 } else { 238 /* Check descent condition */ 239 if (armP->nondescending && *f <= ref - ls->step * fact * ref) break; 240 if (!armP->nondescending && *f <= ref + armP->sigma * gdx) break; 241 ls->step *= armP->beta; 242 } 243 } 244 PetscCall(VecDestroy(&g_old)); 245 246 /* Check termination */ 247 if (PetscIsInfOrNanReal(*f)) { 248 PetscCall(PetscInfo(ls, "Function is inf or nan.\n")); 249 ls->reason = TAOLINESEARCH_FAILED_BADPARAMETER; 250 } else if (ls->step < owlqn_minstep) { 251 PetscCall(PetscInfo(ls, "Step length is below tolerance.\n")); 252 ls->reason = TAOLINESEARCH_HALTED_RTOL; 253 } else if (ls->nfeval >= ls->max_funcs) { 254 PetscCall(PetscInfo(ls, "Number of line search function evals (%" PetscInt_FMT ") > maximum allowed (%" PetscInt_FMT ")\n", ls->nfeval, ls->max_funcs)); 255 ls->reason = TAOLINESEARCH_HALTED_MAXFCN; 256 } 257 if (ls->reason) PetscFunctionReturn(0); 258 259 /* Successful termination, update memory */ 260 ls->reason = TAOLINESEARCH_SUCCESS; 261 armP->lastReference = ref; 262 if (armP->replacementPolicy == REPLACE_FIFO) { 263 armP->memory[armP->current++] = *f; 264 if (armP->current >= armP->memorySize) { armP->current = 0; } 265 } else { 266 armP->current = idx; 267 armP->memory[idx] = *f; 268 } 269 270 /* Update iterate and compute gradient */ 271 PetscCall(VecCopy(armP->work, x)); 272 if (!g_computed) { PetscCall(TaoLineSearchComputeGradient(ls, x, g)); } 273 PetscCall(PetscInfo(ls, "%" PetscInt_FMT " function evals in line search, step = %10.4f\n", ls->nfeval, (double)ls->step)); 274 PetscFunctionReturn(0); 275 } 276 277 /*MC 278 TAOLINESEARCHOWARMIJO - Special line-search type for the Orthant-Wise Limited Quasi-Newton (TAOOWLQN) algorithm. 279 Should not be used with any other algorithm. 280 281 Level: developer 282 283 .keywords: Tao, linesearch 284 M*/ 285 PETSC_EXTERN PetscErrorCode TaoLineSearchCreate_OWArmijo(TaoLineSearch ls) { 286 TaoLineSearch_OWARMIJO *armP; 287 288 PetscFunctionBegin; 289 PetscValidHeaderSpecific(ls, TAOLINESEARCH_CLASSID, 1); 290 PetscCall(PetscNewLog(ls, &armP)); 291 292 armP->memory = NULL; 293 armP->alpha = 1.0; 294 armP->beta = 0.25; 295 armP->beta_inf = 0.25; 296 armP->sigma = 1e-4; 297 armP->memorySize = 1; 298 armP->referencePolicy = REFERENCE_MAX; 299 armP->replacementPolicy = REPLACE_MRU; 300 armP->nondescending = PETSC_FALSE; 301 ls->data = (void *)armP; 302 ls->initstep = 0.1; 303 ls->ops->monitor = NULL; 304 ls->ops->setup = NULL; 305 ls->ops->reset = NULL; 306 ls->ops->apply = TaoLineSearchApply_OWArmijo; 307 ls->ops->view = TaoLineSearchView_OWArmijo; 308 ls->ops->destroy = TaoLineSearchDestroy_OWArmijo; 309 ls->ops->setfromoptions = TaoLineSearchSetFromOptions_OWArmijo; 310 PetscFunctionReturn(0); 311 } 312