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