Actual source code: almmutils.c
1: #include <../src/tao/constrained/impls/almm/almm.h>
2: #include <petsctao.h>
3: #include <petsc/private/petscimpl.h>
4: #include <petsc/private/vecimpl.h>
6: /*@
7: TaoALMMGetType - Retrieve the augmented Lagrangian formulation type for the subproblem.
9: Input Parameter:
10: . tao - the `Tao` context for the `TAOALMM` solver
12: Output Parameter:
13: . type - augmented Lagragrangian type
15: Level: advanced
17: .seealso: `Tao`, `TAOALMM`, `TaoALMMSetType()`, `TaoALMMType`
18: @*/
19: PetscErrorCode TaoALMMGetType(Tao tao, TaoALMMType *type)
20: {
21: PetscFunctionBegin;
23: PetscAssertPointer(type, 2);
24: PetscUseMethod(tao, "TaoALMMGetType_C", (Tao, TaoALMMType *), (tao, type));
25: PetscFunctionReturn(PETSC_SUCCESS);
26: }
28: PetscErrorCode TaoALMMGetType_Private(Tao tao, TaoALMMType *type)
29: {
30: TAO_ALMM *auglag = (TAO_ALMM *)tao->data;
32: PetscFunctionBegin;
33: *type = auglag->type;
34: PetscFunctionReturn(PETSC_SUCCESS);
35: }
37: /*@
38: TaoALMMSetType - Determine the augmented Lagrangian formulation type for the subproblem.
40: Input Parameters:
41: + tao - the `Tao` context for the `TAOALMM` solver
42: - type - augmented Lagragrangian type
44: Level: advanced
46: .seealso: `Tao`, `TAOALMM`, `TaoALMMGetType()`, `TaoALMMType`
47: @*/
48: PetscErrorCode TaoALMMSetType(Tao tao, TaoALMMType type)
49: {
50: PetscFunctionBegin;
52: PetscTryMethod(tao, "TaoALMMSetType_C", (Tao, TaoALMMType), (tao, type));
53: PetscFunctionReturn(PETSC_SUCCESS);
54: }
56: PetscErrorCode TaoALMMSetType_Private(Tao tao, TaoALMMType type)
57: {
58: TAO_ALMM *auglag = (TAO_ALMM *)tao->data;
60: PetscFunctionBegin;
61: PetscCheck(!tao->setupcalled, PetscObjectComm((PetscObject)tao), PETSC_ERR_ORDER, "TaoALMMSetType() must be called before TaoSetUp()");
62: auglag->type = type;
63: PetscFunctionReturn(PETSC_SUCCESS);
64: }
66: /*@
67: TaoALMMGetSubsolver - Retrieve the subsolver being used by `TAOALMM`.
69: Input Parameter:
70: . tao - the `Tao` context for the `TAOALMM` solver
72: Output Parameter:
73: . subsolver - the `Tao` context for the subsolver
75: Level: advanced
77: .seealso: `Tao`, `TAOALMM`, `TaoALMMSetSubsolver()`
78: @*/
79: PetscErrorCode TaoALMMGetSubsolver(Tao tao, Tao *subsolver)
80: {
81: PetscFunctionBegin;
83: PetscAssertPointer(subsolver, 2);
84: PetscUseMethod(tao, "TaoALMMGetSubsolver_C", (Tao, Tao *), (tao, subsolver));
85: PetscFunctionReturn(PETSC_SUCCESS);
86: }
88: PetscErrorCode TaoALMMGetSubsolver_Private(Tao tao, Tao *subsolver)
89: {
90: TAO_ALMM *auglag = (TAO_ALMM *)tao->data;
92: PetscFunctionBegin;
93: *subsolver = auglag->subsolver;
94: PetscFunctionReturn(PETSC_SUCCESS);
95: }
97: /*@
98: TaoALMMSetSubsolver - Changes the subsolver inside `TAOALMM` with the user provided one.
100: Input Parameters:
101: + tao - the `Tao` context for the `TAOALMM` solver
102: - subsolver - the Tao context for the subsolver
104: Level: advanced
106: Note:
107: This is not recommended, instead call `TaoALMMGetSubsolver()` and set the type as desired.
109: .seealso: `Tao`, `TAOALMM`, `TaoALMMGetSubsolver()`
110: @*/
111: PetscErrorCode TaoALMMSetSubsolver(Tao tao, Tao subsolver)
112: {
113: PetscFunctionBegin;
116: PetscTryMethod(tao, "TaoALMMSetSubsolver_C", (Tao, Tao), (tao, subsolver));
117: PetscFunctionReturn(PETSC_SUCCESS);
118: }
120: PetscErrorCode TaoALMMSetSubsolver_Private(Tao tao, Tao subsolver)
121: {
122: TAO_ALMM *auglag = (TAO_ALMM *)tao->data;
123: PetscBool compatible;
125: PetscFunctionBegin;
126: if (subsolver == auglag->subsolver) PetscFunctionReturn(PETSC_SUCCESS);
127: if (tao->bounded) {
128: PetscCall(PetscObjectTypeCompareAny((PetscObject)subsolver, &compatible, TAOSHELL, TAOBNCG, TAOBQNLS, TAOBQNKLS, TAOBQNKTR, TAOBQNKTL, ""));
129: PetscCheck(compatible, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "Subsolver must be a bound-constrained first-order method");
130: } else {
131: PetscCall(PetscObjectTypeCompareAny((PetscObject)subsolver, &compatible, TAOSHELL, TAOCG, TAOLMVM, TAOBNCG, TAOBQNLS, TAOBQNKLS, TAOBQNKTR, TAOBQNKTL, ""));
132: PetscCheck(compatible, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "Subsolver must be a first-order method");
133: }
134: PetscCheck(compatible, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "Subsolver must be a first-order method");
135: PetscCall(PetscObjectReference((PetscObject)subsolver));
136: PetscCall(TaoDestroy(&auglag->subsolver));
137: auglag->subsolver = subsolver;
138: if (tao->setupcalled) {
139: PetscCall(TaoSetSolution(auglag->subsolver, auglag->P));
140: PetscCall(TaoSetObjective(auglag->subsolver, TaoALMMSubsolverObjective_Private, (void *)auglag));
141: PetscCall(TaoSetObjectiveAndGradient(auglag->subsolver, NULL, TaoALMMSubsolverObjectiveAndGradient_Private, (void *)auglag));
142: PetscCall(TaoSetVariableBounds(auglag->subsolver, auglag->PL, auglag->PU));
143: }
144: PetscFunctionReturn(PETSC_SUCCESS);
145: }
147: /*@
148: TaoALMMGetMultipliers - Retrieve a pointer to the Lagrange multipliers.
150: Input Parameter:
151: . tao - the `Tao` context for the `TAOALMM` solver
153: Output Parameter:
154: . Y - vector of Lagrange multipliers
156: Level: advanced
158: Notes:
159: For problems with both equality and inequality constraints,
160: the multipliers are combined together as Y = (Ye, Yi). Users
161: can recover copies of the subcomponents using index sets
162: provided by `TaoALMMGetDualIS()` and use `VecGetSubVector()`.
164: .seealso: `TAOALMM`, `Tao`, `TaoALMMSetMultipliers()`, `TaoALMMGetDualIS()`
165: @*/
166: PetscErrorCode TaoALMMGetMultipliers(Tao tao, Vec *Y)
167: {
168: PetscFunctionBegin;
170: PetscAssertPointer(Y, 2);
171: PetscUseMethod(tao, "TaoALMMGetMultipliers_C", (Tao, Vec *), (tao, Y));
172: PetscFunctionReturn(PETSC_SUCCESS);
173: }
175: PetscErrorCode TaoALMMGetMultipliers_Private(Tao tao, Vec *Y)
176: {
177: TAO_ALMM *auglag = (TAO_ALMM *)tao->data;
179: PetscFunctionBegin;
180: PetscCheck(tao->setupcalled, PetscObjectComm((PetscObject)tao), PETSC_ERR_ORDER, "TaoSetUp() must be called first for scatters to be constructed");
181: *Y = auglag->Y;
182: PetscFunctionReturn(PETSC_SUCCESS);
183: }
185: /*@
186: TaoALMMSetMultipliers - Set user-defined Lagrange multipliers.
188: Input Parameters:
189: + tao - the `Tao` context for the `TAOALMM` solver
190: - Y - vector of Lagrange multipliers
192: Level: advanced
194: Notes:
195: The vector type and parallel layout must match the equality and inequality constraints.
197: The vector must have a local size equal to the sum of the local sizes for the constraint vectors, and a
198: global size equal to the sum of the global sizes of the constraint vectors.
200: This routine is only useful if the user wants to change the
201: parallel distribution of the combined dual vector in problems that
202: feature both equality and inequality constraints. For other tasks,
203: it is strongly recommended that the user retrieve the dual vector
204: created by the solver using `TaoALMMGetMultipliers()`.
206: .seealso: `TAOALMM`, `Tao`, `TaoALMMGetMultipliers()`
207: @*/
208: PetscErrorCode TaoALMMSetMultipliers(Tao tao, Vec Y)
209: {
210: PetscFunctionBegin;
213: PetscTryMethod(tao, "TaoALMMSetMultipliers_C", (Tao, Vec), (tao, Y));
214: PetscFunctionReturn(PETSC_SUCCESS);
215: }
217: PetscErrorCode TaoALMMSetMultipliers_Private(Tao tao, Vec Y)
218: {
219: TAO_ALMM *auglag = (TAO_ALMM *)tao->data;
220: VecType Ytype;
221: PetscInt Nuser, Neq, Nineq, N;
222: PetscBool same = PETSC_FALSE;
224: PetscFunctionBegin;
225: /* no-op if user provides vector from TaoALMMGetMultipliers() */
226: if (Y == auglag->Y) PetscFunctionReturn(PETSC_SUCCESS);
227: /* make sure vector type is same as equality and inequality constraints */
228: if (tao->eq_constrained) {
229: PetscCall(VecGetType(tao->constraints_equality, &Ytype));
230: } else {
231: PetscCall(VecGetType(tao->constraints_inequality, &Ytype));
232: }
233: PetscCall(PetscObjectTypeCompare((PetscObject)Y, Ytype, &same));
234: PetscCheck(same, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "Given vector for multipliers is not the same type as constraint vectors");
235: /* make sure global size matches sum of equality and inequality */
236: if (tao->eq_constrained) {
237: PetscCall(VecGetSize(tao->constraints_equality, &Neq));
238: } else {
239: Neq = 0;
240: }
241: if (tao->ineq_constrained) {
242: PetscCall(VecGetSize(tao->constraints_inequality, &Nineq));
243: } else {
244: Nineq = 0;
245: }
246: N = Neq + Nineq;
247: PetscCall(VecGetSize(Y, &Nuser));
248: PetscCheck(Nuser == N, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "Given vector has wrong global size");
249: /* if there is only one type of constraint, then we need the local size to match too */
250: if (Neq == 0) {
251: PetscCall(VecGetLocalSize(tao->constraints_inequality, &Nineq));
252: PetscCall(VecGetLocalSize(Y, &Nuser));
253: PetscCheck(Nuser == Nineq, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "Given vector has wrong local size");
254: }
255: if (Nineq == 0) {
256: PetscCall(VecGetLocalSize(tao->constraints_equality, &Neq));
257: PetscCall(VecGetLocalSize(Y, &Nuser));
258: PetscCheck(Nuser == Neq, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_INCOMP, "Given vector has wrong local size");
259: }
260: /* if we got here, the given vector is compatible so we can replace the current one */
261: PetscCall(PetscObjectReference((PetscObject)Y));
262: PetscCall(VecDestroy(&auglag->Y));
263: auglag->Y = Y;
264: /* if there are both types of constraints and the solver has already been set up,
265: then we need to regenerate VecScatter objects for the new combined dual vector */
266: if (tao->setupcalled && tao->eq_constrained && tao->ineq_constrained) {
267: PetscCall(VecDestroy(&auglag->C));
268: PetscCall(VecDuplicate(auglag->Y, &auglag->C));
269: PetscCall(VecScatterDestroy(&auglag->Yscatter[0]));
270: PetscCall(VecScatterCreate(auglag->Y, auglag->Yis[0], auglag->Ye, NULL, &auglag->Yscatter[0]));
271: PetscCall(VecScatterDestroy(&auglag->Yscatter[1]));
272: PetscCall(VecScatterCreate(auglag->Y, auglag->Yis[1], auglag->Yi, NULL, &auglag->Yscatter[1]));
273: }
274: PetscFunctionReturn(PETSC_SUCCESS);
275: }
277: /*@
278: TaoALMMGetPrimalIS - Retrieve the index set that identifies optimization
279: and slack variable components of the subsolver's solution vector.
281: Input Parameter:
282: . tao - the `Tao` context for the `TAOALMM` solver
284: Output Parameters:
285: + opt_is - index set associated with the optimization variables (`NULL` if not needed)
286: - slack_is - index set associated with the slack variables (`NULL` if not needed)
288: Level: advanced
290: .seealso: `TAOALMM`, `Tao`, `IS`, `TaoALMMGetPrimalVector()`
291: @*/
292: PetscErrorCode TaoALMMGetPrimalIS(Tao tao, IS *opt_is, IS *slack_is)
293: {
294: PetscFunctionBegin;
296: PetscUseMethod(tao, "TaoALMMGetPrimalIS_C", (Tao, IS *, IS *), (tao, opt_is, slack_is));
297: PetscFunctionReturn(PETSC_SUCCESS);
298: }
300: PetscErrorCode TaoALMMGetPrimalIS_Private(Tao tao, IS *opt_is, IS *slack_is)
301: {
302: TAO_ALMM *auglag = (TAO_ALMM *)tao->data;
304: PetscFunctionBegin;
305: PetscCheck(tao->ineq_constrained, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_WRONGSTATE, "Primal space has index sets only for inequality constrained problems");
306: PetscCheck(tao->setupcalled, PetscObjectComm((PetscObject)tao), PETSC_ERR_ORDER, "TaoSetUp() must be called first for index sets to be constructed");
307: if (opt_is) *opt_is = auglag->Pis[0];
308: if (slack_is) *slack_is = auglag->Pis[1];
309: PetscFunctionReturn(PETSC_SUCCESS);
310: }
312: /*@
313: TaoALMMGetDualIS - Retrieve the index set that identifies equality
314: and inequality constraint components of the dual vector returned
315: by `TaoALMMGetMultipliers()`.
317: Input Parameter:
318: . tao - the Tao context for the `TAOALMM` solver
320: Output Parameters:
321: + eq_is - index set associated with the equality constraints (`NULL` if not needed)
322: - ineq_is - index set associated with the inequality constraints (`NULL` if not needed)
324: Level: advanced
326: .seealso: `TAOALMM`, `Tao`, `TaoALMMGetMultipliers()`
327: @*/
328: PetscErrorCode TaoALMMGetDualIS(Tao tao, IS *eq_is, IS *ineq_is)
329: {
330: PetscFunctionBegin;
332: PetscUseMethod(tao, "TaoALMMGetDualIS_C", (Tao, IS *, IS *), (tao, eq_is, ineq_is));
333: PetscFunctionReturn(PETSC_SUCCESS);
334: }
336: PetscErrorCode TaoALMMGetDualIS_Private(Tao tao, IS *eq_is, IS *ineq_is)
337: {
338: TAO_ALMM *auglag = (TAO_ALMM *)tao->data;
340: PetscFunctionBegin;
342: PetscCheck(tao->ineq_constrained && tao->ineq_constrained, PetscObjectComm((PetscObject)tao), PETSC_ERR_ARG_WRONGSTATE, "Dual space has index sets only when problem has both equality and inequality constraints");
343: PetscCheck(tao->setupcalled, PetscObjectComm((PetscObject)tao), PETSC_ERR_ORDER, "TaoSetUp() must be called first for index sets to be constructed");
344: if (eq_is) *eq_is = auglag->Yis[0];
345: if (ineq_is) *ineq_is = auglag->Yis[1];
346: PetscFunctionReturn(PETSC_SUCCESS);
347: }