Actual source code: taosolver_bounds.c

  1: #include <petsc/private/taoimpl.h>

  3: /*@
  4:   TaoSetVariableBounds - Sets the upper and lower bounds

  6:   Logically collective on Tao

  8:   Input Parameters:
  9: + tao - the Tao context
 10: . XL  - vector of lower bounds
 11: - XU  - vector of upper bounds

 13:   Level: beginner

 15: .seealso: TaoSetObjectiveRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine()
 16: @*/

 18: PetscErrorCode TaoSetVariableBounds(Tao tao, Vec XL, Vec XU)
 19: {

 24:   if (XL) {
 26:     PetscObjectReference((PetscObject)XL);
 27:   }
 28:   if (XU) {
 30:     PetscObjectReference((PetscObject)XU);
 31:   }
 32:   VecDestroy(&tao->XL);
 33:   VecDestroy(&tao->XU);
 34:   tao->XL = XL;
 35:   tao->XU = XU;
 36:   tao->bounded = PETSC_TRUE;
 37:   return(0);
 38: }

 40: /*@C
 41:   TaoSetVariableBoundsRoutine - Sets a function to be used to compute variable bounds

 43:   Logically collective on Tao

 45:   Input Parameters:
 46: + tao - the Tao context
 47: . func - the bounds computation routine
 48: - ctx - [optional] user-defined context for private data for the bounds computation (may be NULL)

 50:   Calling sequence of func:
 51: $      func (Tao tao, Vec xl, Vec xu);

 53: + tao - the Tao
 54: . xl  - vector of lower bounds
 55: . xu  - vector of upper bounds
 56: - ctx - the (optional) user-defined function context

 58:   Level: beginner

 60: .seealso: TaoSetObjectiveRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine(), TaoSetVariableBounds()

 62: Note: The func passed in to TaoSetVariableBoundsRoutine() takes
 63: precedence over any values set in TaoSetVariableBounds().

 65: @*/
 66: PetscErrorCode TaoSetVariableBoundsRoutine(Tao tao, PetscErrorCode (*func)(Tao, Vec, Vec, void*), void *ctx)
 67: {
 70:   tao->user_boundsP = ctx;
 71:   tao->ops->computebounds = func;
 72:   tao->bounded = PETSC_TRUE;
 73:   return(0);
 74: }

 76: PetscErrorCode TaoGetVariableBounds(Tao tao, Vec *XL, Vec *XU)
 77: {
 80:   if (XL) {
 81:     *XL=tao->XL;
 82:   }
 83:   if (XU) {
 84:     *XU=tao->XU;
 85:   }
 86:   return(0);
 87: }

 89: /*@C
 90:    TaoComputeVariableBounds - Compute the variable bounds using the
 91:    routine set by TaoSetVariableBoundsRoutine().

 93:    Collective on Tao

 95:    Input Parameters:
 96: .  tao - the Tao context

 98:    Level: developer

100: .seealso: TaoSetVariableBoundsRoutine(), TaoSetVariableBounds()
101: @*/

103: PetscErrorCode TaoComputeVariableBounds(Tao tao)
104: {

109:   PetscStackPush("Tao compute variable bounds");
110:   if (!tao->XL || !tao->XU) {
111:     if (!tao->solution) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ORDER,"TaoSetInitialVector must be called before TaoComputeVariableBounds");
112:     VecDuplicate(tao->solution, &tao->XL);
113:     VecSet(tao->XL, PETSC_NINFINITY);
114:     VecDuplicate(tao->solution, &tao->XU);
115:     VecSet(tao->XU, PETSC_INFINITY);
116:   }
117:   if (tao->ops->computebounds) {
118:     (*tao->ops->computebounds)(tao,tao->XL,tao->XU,tao->user_boundsP);
119:   }
120:   PetscStackPop;
121:   return(0);
122: }

124: /*@
125:   TaoSetInequalityBounds - Sets the upper and lower bounds

127:   Logically collective on Tao

129:   Input Parameters:
130: + tao - the Tao context
131: . IL  - vector of lower bounds
132: - IU  - vector of upper bounds

134:   Level: beginner

136: .seealso: TaoSetObjectiveRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine()
137: @*/

139: PetscErrorCode TaoSetInequalityBounds(Tao tao, Vec IL, Vec IU)
140: {

145:   if (IL) {
147:     PetscObjectReference((PetscObject)IL);
148:   }
149:   if (IU) {
151:     PetscObjectReference((PetscObject)IU);
152:   }
153:   VecDestroy(&tao->IL);
154:   VecDestroy(&tao->IU);
155:   tao->IL = IL;
156:   tao->IU = IU;
157:   tao->ineq_doublesided = PETSC_TRUE;
158:   return(0);
159: }


162: PetscErrorCode TaoGetInequalityBounds(Tao tao, Vec *IL, Vec *IU)
163: {
166:   if (IL) {
167:     *IL=tao->IL;
168:   }
169:   if (IU) {
170:     *IU=tao->IU;
171:   }
172:   return(0);
173: }

175: /*@C
176:    TaoComputeConstraints - Compute the variable bounds using the
177:    routine set by TaoSetConstraintsRoutine().

179:    Collective on Tao

181:    Input Parameters:
182: .  tao - the Tao context

184:    Level: developer

186: .seealso: TaoSetConstraintsRoutine(), TaoComputeJacobian()
187: @*/

189: PetscErrorCode TaoComputeConstraints(Tao tao, Vec X, Vec C)
190: {


200:   if (!tao->ops->computeconstraints) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"TaoSetConstraintsRoutine() has not been called");
201:   if (!tao->solution) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"TaoSetInitialVector must be called before TaoComputeConstraints");
202:   PetscLogEventBegin(TAO_ConstraintsEval,tao,X,C,NULL);
203:   PetscStackPush("Tao constraints evaluation routine");
204:   (*tao->ops->computeconstraints)(tao,X,C,tao->user_conP);
205:   PetscStackPop;
206:   PetscLogEventEnd(TAO_ConstraintsEval,tao,X,C,NULL);
207:   tao->nconstraints++;
208:   return(0);
209: }

211: /*@C
212:   TaoSetConstraintsRoutine - Sets a function to be used to compute constraints.  TAO only handles constraints under certain conditions, see manual for details

214:   Logically collective on Tao

216:   Input Parameters:
217: + tao - the Tao context
218: . c   - A vector that will be used to store constraint evaluation
219: . func - the bounds computation routine
220: - ctx - [optional] user-defined context for private data for the constraints computation (may be NULL)

222:   Calling sequence of func:
223: $      func (Tao tao, Vec x, Vec c, void *ctx);

225: + tao - the Tao
226: . x   - point to evaluate constraints
227: . c   - vector constraints evaluated at x
228: - ctx - the (optional) user-defined function context

230:   Level: intermediate

232: .seealso: TaoSetObjectiveRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine(), TaoSetVariablevBounds()

234: @*/
235: PetscErrorCode TaoSetConstraintsRoutine(Tao tao, Vec c, PetscErrorCode (*func)(Tao, Vec, Vec, void*), void *ctx)
236: {
239:     tao->constrained = PETSC_TRUE;
240:     tao->constraints = c;
241:     tao->user_conP = ctx;
242:     tao->ops->computeconstraints = func;
243:     return(0);
244: }

246: /*@
247:   TaoComputeDualVariables - Computes the dual vectors corresponding to the bounds
248:   of the variables

250:   Collective on Tao

252:   Input Parameters:
253: . tao - the Tao context

255:   Output Parameter:
256: + DL - dual variable vector for the lower bounds
257: - DU - dual variable vector for the upper bounds

259:   Level: advanced

261:   Note:
262:   DL and DU should be created before calling this routine.  If calling
263:   this routine after using an unconstrained solver, DL and DU are set to all
264:   zeros.

266:   Level: advanced

268: .seealso: TaoComputeObjective(), TaoSetVariableBounds()
269: @*/
270: PetscErrorCode TaoComputeDualVariables(Tao tao, Vec DL, Vec DU)
271: {
279:   if (tao->ops->computedual) {
280:     (*tao->ops->computedual)(tao,DL,DU);
281:   }  else {
282:     VecSet(DL,0.0);
283:     VecSet(DU,0.0);
284:   }
285:   return(0);
286: }

288: /*@
289:   TaoGetDualVariables - Gets pointers to the dual vectors

291:   Collective on Tao

293:   Input Parameters:
294: . tao - the Tao context

296:   Output Parameter:
297: + DE - dual variable vector for the lower bounds
298: - DI - dual variable vector for the upper bounds

300:   Level: advanced

302: .seealso: TaoComputeDualVariables()
303: @*/
304: PetscErrorCode TaoGetDualVariables(Tao tao, Vec *DE, Vec *DI)
305: {
308:   if (DE) {
309:     *DE = tao->DE;
310:   }
311:   if (DI) {
312:     *DI = tao->DI;
313:   }
314:   return(0);
315: }

317: /*@C
318:   TaoSetEqualityConstraintsRoutine - Sets a function to be used to compute constraints.  TAO only handles constraints under certain conditions, see manual for details

320:   Logically collective on Tao

322:   Input Parameters:
323: + tao - the Tao context
324: . ce   - A vector that will be used to store equality constraint evaluation
325: . func - the bounds computation routine
326: - ctx - [optional] user-defined context for private data for the equality constraints computation (may be NULL)

328:   Calling sequence of func:
329: $      func (Tao tao, Vec x, Vec ce, void *ctx);

331: + tao - the Tao
332: . x   - point to evaluate equality constraints
333: . ce   - vector of equality constraints evaluated at x
334: - ctx - the (optional) user-defined function context

336:   Level: intermediate

338: .seealso: TaoSetObjectiveRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine(), TaoSetVariableBounds()

340: @*/
341: PetscErrorCode TaoSetEqualityConstraintsRoutine(Tao tao, Vec ce, PetscErrorCode (*func)(Tao, Vec, Vec, void*), void *ctx)
342: {

347:   if (ce) {
349:     PetscObjectReference((PetscObject)ce);
350:   }
351:   VecDestroy(&tao->constraints_equality);
352:   tao->eq_constrained = PETSC_TRUE;
353:   tao->constraints_equality = ce;
354:   tao->user_con_equalityP = ctx;
355:   tao->ops->computeequalityconstraints = func;
356:   return(0);
357: }


360: /*@C
361:   TaoSetInequalityConstraintsRoutine - Sets a function to be used to compute constraints.  TAO only handles constraints under certain conditions, see manual for details

363:   Logically collective on Tao

365:   Input Parameters:
366: + tao - the Tao context
367: . ci   - A vector that will be used to store inequality constraint evaluation
368: . func - the bounds computation routine
369: - ctx - [optional] user-defined context for private data for the inequality constraints computation (may be NULL)

371:   Calling sequence of func:
372: $      func (Tao tao, Vec x, Vec ci, void *ctx);

374: + tao - the Tao
375: . x   - point to evaluate inequality constraints
376: . ci   - vector of inequality constraints evaluated at x
377: - ctx - the (optional) user-defined function context

379:   Level: intermediate

381: .seealso: TaoSetObjectiveRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine(), TaoSetVariableBounds()

383: @*/
384: PetscErrorCode TaoSetInequalityConstraintsRoutine(Tao tao, Vec ci, PetscErrorCode (*func)(Tao, Vec, Vec, void*), void *ctx)
385: {

390:   if (ci) {
392:     PetscObjectReference((PetscObject)ci);
393:   }
394:   VecDestroy(&tao->constraints_inequality);
395:   tao->constraints_inequality = ci;
396:   tao->ineq_constrained = PETSC_TRUE;
397:   tao->user_con_inequalityP = ctx;
398:   tao->ops->computeinequalityconstraints = func;
399:   return(0);
400: }


403: /*@C
404:    TaoComputeEqualityConstraints - Compute the variable bounds using the
405:    routine set by TaoSetEqualityConstraintsRoutine().

407:    Collective on Tao

409:    Input Parameters:
410: .  tao - the Tao context

412:    Level: developer

414: .seealso: TaoSetEqualityConstraintsRoutine(), TaoComputeJacobianEquality()
415: @*/

417: PetscErrorCode TaoComputeEqualityConstraints(Tao tao, Vec X, Vec CE)
418: {


428:   if (!tao->ops->computeequalityconstraints) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"TaoSetEqualityConstraintsRoutine() has not been called");
429:   if (!tao->solution) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"TaoSetInitialVector must be called before TaoComputeEqualityConstraints");
430:   PetscLogEventBegin(TAO_ConstraintsEval,tao,X,CE,NULL);
431:   PetscStackPush("Tao equality constraints evaluation routine");
432:   (*tao->ops->computeequalityconstraints)(tao,X,CE,tao->user_con_equalityP);
433:   PetscStackPop;
434:   PetscLogEventEnd(TAO_ConstraintsEval,tao,X,CE,NULL);
435:   tao->nconstraints++;
436:   return(0);
437: }


440: /*@C
441:    TaoComputeInequalityConstraints - Compute the variable bounds using the
442:    routine set by TaoSetInequalityConstraintsRoutine().

444:    Collective on Tao

446:    Input Parameters:
447: .  tao - the Tao context

449:    Level: developer

451: .seealso: TaoSetInequalityConstraintsRoutine(), TaoComputeJacobianInequality()
452: @*/

454: PetscErrorCode TaoComputeInequalityConstraints(Tao tao, Vec X, Vec CI)
455: {


465:   if (!tao->ops->computeinequalityconstraints) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"TaoSetInequalityConstraintsRoutine() has not been called");
466:   if (!tao->solution) SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"TaoSetInitialVector must be called before TaoComputeInequalityConstraints");
467:   PetscLogEventBegin(TAO_ConstraintsEval,tao,X,CI,NULL);
468:   PetscStackPush("Tao inequality constraints evaluation routine");
469:   (*tao->ops->computeinequalityconstraints)(tao,X,CI,tao->user_con_inequalityP);
470:   PetscStackPop;
471:   PetscLogEventEnd(TAO_ConstraintsEval,tao,X,CI,NULL);
472:   tao->nconstraints++;
473:   return(0);
474: }