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