Actual source code: ex1.c
1: /* Program usage: mpiexec -n 2 ./ex1 [-help] [all TAO options] */
3: /* ----------------------------------------------------------------------
4: min f(x) = (x0 - 2)^2 + (x1 - 2)^2 - 2*(x0 + x1)
5: s.t. x0^2 + x1 - 2 = 0
6: 0 <= x0^2 - x1 <= 1
7: -1 <= x0, x1 <= 2
8: -->
9: g(x) = 0
10: h(x) >= 0
11: -1 <= x0, x1 <= 2
12: where
13: g(x) = x0^2 + x1 - 2
14: h(x) = [x0^2 - x1
15: 1 -(x0^2 - x1)]
16: ---------------------------------------------------------------------- */
18: #include <petsctao.h>
20: static char help[] = "Solves constrained optimization problem using pdipm.\n\
21: Input parameters include:\n\
22: -tao_type pdipm : sets Tao solver\n\
23: -no_eq : removes the equaility constraints from the problem\n\
24: -init_view : view initial object setup\n\
25: -snes_fd : snes with finite difference Jacobian (needed for pdipm)\n\
26: -snes_compare_explicit : compare user Jacobian with finite difference Jacobian \n\
27: -tao_monitor_constraint_norm : convergence monitor with constraint norm \n\
28: -tao_monitor_solution : view exact solution at each iteration\n\
29: Note: external package MUMPS is required to run pdipm in parallel. This is designed for a maximum of 2 processors, the code will error if size > 2.\n";
31: /*
32: User-defined application context - contains data needed by the application
33: */
34: typedef struct {
35: PetscInt n; /* Global length of x */
36: PetscInt ne; /* Global number of equality constraints */
37: PetscInt ni; /* Global number of inequality constraints */
38: PetscBool noeqflag, initview;
39: Vec x, xl, xu;
40: Vec ce, ci, bl, bu, Xseq;
41: Mat Ae, Ai, H;
42: VecScatter scat;
43: } AppCtx;
45: /* -------- User-defined Routines --------- */
46: PetscErrorCode InitializeProblem(AppCtx *);
47: PetscErrorCode DestroyProblem(AppCtx *);
48: PetscErrorCode FormFunctionGradient(Tao, Vec, PetscReal *, Vec, void *);
49: PetscErrorCode FormPDIPMHessian(Tao, Vec, Mat, Mat, void *);
50: PetscErrorCode FormInequalityConstraints(Tao, Vec, Vec, void *);
51: PetscErrorCode FormEqualityConstraints(Tao, Vec, Vec, void *);
52: PetscErrorCode FormInequalityJacobian(Tao, Vec, Mat, Mat, void *);
53: PetscErrorCode FormEqualityJacobian(Tao, Vec, Mat, Mat, void *);
55: int main(int argc, char **argv)
56: {
57: Tao tao;
58: KSP ksp;
59: PC pc;
60: AppCtx user; /* application context */
61: Vec x, G, CI, CE;
62: PetscMPIInt size;
63: TaoType type;
64: PetscReal f;
65: PetscBool pdipm, mumps;
67: PetscFunctionBeginUser;
68: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
69: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
70: PetscCheck(size <= 2, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "More than 2 processors detected. Example written to use max of 2 processors.");
72: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "---- Constrained Problem -----\n"));
73: PetscCall(InitializeProblem(&user)); /* sets up problem, function below */
75: PetscCall(TaoCreate(PETSC_COMM_WORLD, &tao));
76: PetscCall(TaoSetType(tao, TAOALMM));
77: PetscCall(TaoSetSolution(tao, user.x));
78: PetscCall(TaoSetVariableBounds(tao, user.xl, user.xu));
79: PetscCall(TaoSetObjectiveAndGradient(tao, NULL, FormFunctionGradient, (void *)&user));
80: PetscCall(TaoSetTolerances(tao, 1.e-4, 0.0, 0.0));
81: PetscCall(TaoSetConstraintTolerances(tao, 1.e-4, 0.0));
82: PetscCall(TaoSetMaximumFunctionEvaluations(tao, 1e6));
83: PetscCall(TaoSetFromOptions(tao));
85: if (!user.noeqflag) {
86: PetscCall(TaoSetEqualityConstraintsRoutine(tao, user.ce, FormEqualityConstraints, (void *)&user));
87: PetscCall(TaoSetJacobianEqualityRoutine(tao, user.Ae, user.Ae, FormEqualityJacobian, (void *)&user));
88: }
89: PetscCall(TaoSetInequalityConstraintsRoutine(tao, user.ci, FormInequalityConstraints, (void *)&user));
90: PetscCall(TaoSetJacobianInequalityRoutine(tao, user.Ai, user.Ai, FormInequalityJacobian, (void *)&user));
92: PetscCall(TaoGetType(tao, &type));
93: PetscCall(PetscObjectTypeCompare((PetscObject)tao, TAOPDIPM, &pdipm));
94: if (pdipm) {
95: /*
96: PDIPM produces an indefinite KKT matrix with some zeros along the diagonal
97: Inverting this indefinite matrix requires PETSc to be configured with MUMPS
98: */
99: PetscCall(PetscHasExternalPackage("mumps", &mumps));
100: PetscCheck(mumps, PetscObjectComm((PetscObject)tao), PETSC_ERR_SUP, "TAOPDIPM requires PETSc to be configured with MUMPS (--download-mumps)");
101: PetscCall(TaoGetKSP(tao, &ksp));
102: PetscCall(KSPGetPC(ksp, &pc));
103: PetscCall(PCFactorSetMatSolverType(pc, MATSOLVERMUMPS));
104: PetscCall(KSPSetType(ksp, KSPPREONLY));
105: PetscCall(PCSetType(pc, PCCHOLESKY));
106: PetscCall(KSPSetFromOptions(ksp));
107: PetscCall(TaoSetHessian(tao, user.H, user.H, FormPDIPMHessian, (void *)&user));
108: }
110: /* Print out an initial view of the problem */
111: if (user.initview) {
112: PetscCall(TaoSetUp(tao));
113: PetscCall(VecDuplicate(user.x, &G));
114: PetscCall(FormFunctionGradient(tao, user.x, &f, G, (void *)&user));
115: PetscCall(FormPDIPMHessian(tao, user.x, user.H, user.H, (void *)&user));
116: PetscCall(PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_WORLD));
117: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial point X:\n"));
118: PetscCall(VecView(user.x, PETSC_VIEWER_STDOUT_WORLD));
119: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial objective f(x) = %g\n", (double)f));
120: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial gradient and Hessian:\n"));
121: PetscCall(VecView(G, PETSC_VIEWER_STDOUT_WORLD));
122: PetscCall(MatView(user.H, PETSC_VIEWER_STDOUT_WORLD));
123: PetscCall(VecDestroy(&G));
124: PetscCall(FormInequalityJacobian(tao, user.x, user.Ai, user.Ai, (void *)&user));
125: PetscCall(MatCreateVecs(user.Ai, NULL, &CI));
126: PetscCall(FormInequalityConstraints(tao, user.x, CI, (void *)&user));
127: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial inequality constraints and Jacobian:\n"));
128: PetscCall(VecView(CI, PETSC_VIEWER_STDOUT_WORLD));
129: PetscCall(MatView(user.Ai, PETSC_VIEWER_STDOUT_WORLD));
130: PetscCall(VecDestroy(&CI));
131: if (!user.noeqflag) {
132: PetscCall(FormEqualityJacobian(tao, user.x, user.Ae, user.Ae, (void *)&user));
133: PetscCall(MatCreateVecs(user.Ae, NULL, &CE));
134: PetscCall(FormEqualityConstraints(tao, user.x, CE, (void *)&user));
135: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial equality constraints and Jacobian:\n"));
136: PetscCall(VecView(CE, PETSC_VIEWER_STDOUT_WORLD));
137: PetscCall(MatView(user.Ae, PETSC_VIEWER_STDOUT_WORLD));
138: PetscCall(VecDestroy(&CE));
139: }
140: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\n"));
141: PetscCall(PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_WORLD));
142: }
144: PetscCall(TaoSolve(tao));
145: PetscCall(TaoGetSolution(tao, &x));
146: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Found solution:\n"));
147: PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD));
149: /* Free objects */
150: PetscCall(DestroyProblem(&user));
151: PetscCall(TaoDestroy(&tao));
152: PetscCall(PetscFinalize());
153: return PETSC_SUCCESS;
154: }
156: PetscErrorCode InitializeProblem(AppCtx *user)
157: {
158: PetscMPIInt size;
159: PetscMPIInt rank;
160: PetscInt nloc, neloc, niloc;
162: PetscFunctionBegin;
163: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
164: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
165: user->noeqflag = PETSC_FALSE;
166: PetscCall(PetscOptionsGetBool(NULL, NULL, "-no_eq", &user->noeqflag, NULL));
167: user->initview = PETSC_FALSE;
168: PetscCall(PetscOptionsGetBool(NULL, NULL, "-init_view", &user->initview, NULL));
170: /* Tell user the correct solution, not an error checking */
171: if (!user->noeqflag) {
172: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Expected solution: f(1, 1) = -2\n"));
173: } else {
174: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Expected solution (-no_eq): f(1.73205, 2) = -7.3923\n"));
175: }
177: /* create vector x and set initial values */
178: user->n = 2; /* global length */
179: nloc = (size == 1) ? user->n : 1;
180: PetscCall(VecCreate(PETSC_COMM_WORLD, &user->x));
181: PetscCall(VecSetSizes(user->x, nloc, user->n));
182: PetscCall(VecSetFromOptions(user->x));
183: PetscCall(VecSet(user->x, 0.0));
185: /* create and set lower and upper bound vectors */
186: PetscCall(VecDuplicate(user->x, &user->xl));
187: PetscCall(VecDuplicate(user->x, &user->xu));
188: PetscCall(VecSet(user->xl, -1.0));
189: PetscCall(VecSet(user->xu, 2.0));
191: /* create scater to zero */
192: PetscCall(VecScatterCreateToZero(user->x, &user->scat, &user->Xseq));
194: user->ne = 1;
195: user->ni = 2;
196: neloc = (rank == 0) ? user->ne : 0;
197: niloc = (size == 1) ? user->ni : 1;
199: if (!user->noeqflag) {
200: PetscCall(VecCreate(PETSC_COMM_WORLD, &user->ce)); /* a 1x1 vec for equality constraints */
201: PetscCall(VecSetSizes(user->ce, neloc, user->ne));
202: PetscCall(VecSetFromOptions(user->ce));
203: PetscCall(VecSetUp(user->ce));
204: }
206: PetscCall(VecCreate(PETSC_COMM_WORLD, &user->ci)); /* a 2x1 vec for inequality constraints */
207: PetscCall(VecSetSizes(user->ci, niloc, user->ni));
208: PetscCall(VecSetFromOptions(user->ci));
209: PetscCall(VecSetUp(user->ci));
211: /* nexn & nixn matrices for equally and inequalty constraints */
212: if (!user->noeqflag) {
213: PetscCall(MatCreate(PETSC_COMM_WORLD, &user->Ae));
214: PetscCall(MatSetSizes(user->Ae, neloc, nloc, user->ne, user->n));
215: PetscCall(MatSetFromOptions(user->Ae));
216: PetscCall(MatSetUp(user->Ae));
217: }
219: PetscCall(MatCreate(PETSC_COMM_WORLD, &user->Ai));
220: PetscCall(MatSetSizes(user->Ai, niloc, nloc, user->ni, user->n));
221: PetscCall(MatSetFromOptions(user->Ai));
222: PetscCall(MatSetUp(user->Ai));
224: PetscCall(MatCreate(PETSC_COMM_WORLD, &user->H));
225: PetscCall(MatSetSizes(user->H, nloc, nloc, user->n, user->n));
226: PetscCall(MatSetFromOptions(user->H));
227: PetscCall(MatSetUp(user->H));
228: PetscFunctionReturn(PETSC_SUCCESS);
229: }
231: PetscErrorCode DestroyProblem(AppCtx *user)
232: {
233: PetscFunctionBegin;
234: if (!user->noeqflag) PetscCall(MatDestroy(&user->Ae));
235: PetscCall(MatDestroy(&user->Ai));
236: PetscCall(MatDestroy(&user->H));
238: PetscCall(VecDestroy(&user->x));
239: if (!user->noeqflag) PetscCall(VecDestroy(&user->ce));
240: PetscCall(VecDestroy(&user->ci));
241: PetscCall(VecDestroy(&user->xl));
242: PetscCall(VecDestroy(&user->xu));
243: PetscCall(VecDestroy(&user->Xseq));
244: PetscCall(VecScatterDestroy(&user->scat));
245: PetscFunctionReturn(PETSC_SUCCESS);
246: }
248: /* Evaluate
249: f(x) = (x0 - 2)^2 + (x1 - 2)^2 - 2*(x0 + x1)
250: G = grad f = [2*(x0 - 2) - 2, 2*(x1 - 2) - 2] = 2*X - 6
251: */
252: PetscErrorCode FormFunctionGradient(Tao tao, Vec X, PetscReal *f, Vec G, void *ctx)
253: {
254: const PetscScalar *x;
255: MPI_Comm comm;
256: PetscMPIInt rank;
257: PetscReal fin;
258: AppCtx *user = (AppCtx *)ctx;
259: Vec Xseq = user->Xseq;
260: VecScatter scat = user->scat;
262: PetscFunctionBegin;
263: /* f = (x0 - 2)^2 + (x1 - 2)^2 - 2*(x0 + x1) */
264: PetscCall(PetscObjectGetComm((PetscObject)tao, &comm));
265: PetscCallMPI(MPI_Comm_rank(comm, &rank));
267: PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD));
268: PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD));
270: fin = 0.0;
271: if (rank == 0) {
272: PetscCall(VecGetArrayRead(Xseq, &x));
273: fin = (x[0] - 2.0) * (x[0] - 2.0) + (x[1] - 2.0) * (x[1] - 2.0) - 2.0 * (x[0] + x[1]);
274: PetscCall(VecRestoreArrayRead(Xseq, &x));
275: }
276: PetscCall(MPIU_Allreduce(&fin, f, 1, MPIU_REAL, MPIU_SUM, comm));
278: /* G = 2*X - 6 */
279: PetscCall(VecSet(G, -6.0));
280: PetscCall(VecAXPY(G, 2.0, X));
281: PetscFunctionReturn(PETSC_SUCCESS);
282: }
284: /* Evaluate PDIPM Hessian, see Eqn(22) in http://doi.org/10.1049/gtd2.12708
285: H = Wxx = fxx + Jacobian(grad g^T*DE) - Jacobian(grad h^T*DI)]
286: = fxx + Jacobin([2*x0; 1]de[0]) - Jacobian([2*x0, -2*x0; -1, 1][di[0] di[1]]^T)
287: = [2 0; 0 2] + [2*de[0] 0; 0 0] - [2*di[0]-2*di[1], 0; 0, 0]
288: = [ 2*(1+de[0]-di[0]+di[1]), 0;
289: 0, 2]
290: */
291: PetscErrorCode FormPDIPMHessian(Tao tao, Vec x, Mat H, Mat Hpre, void *ctx)
292: {
293: AppCtx *user = (AppCtx *)ctx;
294: Vec DE, DI;
295: const PetscScalar *de, *di;
296: PetscInt zero = 0, one = 1;
297: PetscScalar two = 2.0;
298: PetscScalar val = 0.0;
299: Vec Deseq, Diseq;
300: VecScatter Descat, Discat;
301: PetscMPIInt rank;
302: MPI_Comm comm;
304: PetscFunctionBegin;
305: PetscCall(TaoGetDualVariables(tao, &DE, &DI));
307: PetscCall(PetscObjectGetComm((PetscObject)tao, &comm));
308: PetscCallMPI(MPI_Comm_rank(comm, &rank));
310: if (!user->noeqflag) {
311: PetscCall(VecScatterCreateToZero(DE, &Descat, &Deseq));
312: PetscCall(VecScatterBegin(Descat, DE, Deseq, INSERT_VALUES, SCATTER_FORWARD));
313: PetscCall(VecScatterEnd(Descat, DE, Deseq, INSERT_VALUES, SCATTER_FORWARD));
314: }
315: PetscCall(VecScatterCreateToZero(DI, &Discat, &Diseq));
316: PetscCall(VecScatterBegin(Discat, DI, Diseq, INSERT_VALUES, SCATTER_FORWARD));
317: PetscCall(VecScatterEnd(Discat, DI, Diseq, INSERT_VALUES, SCATTER_FORWARD));
319: if (rank == 0) {
320: if (!user->noeqflag) { PetscCall(VecGetArrayRead(Deseq, &de)); /* places equality constraint dual into array */ }
321: PetscCall(VecGetArrayRead(Diseq, &di)); /* places inequality constraint dual into array */
323: if (!user->noeqflag) {
324: val = 2.0 * (1 + de[0] - di[0] + di[1]);
325: PetscCall(VecRestoreArrayRead(Deseq, &de));
326: PetscCall(VecRestoreArrayRead(Diseq, &di));
327: } else {
328: val = 2.0 * (1 - di[0] + di[1]);
329: }
330: PetscCall(VecRestoreArrayRead(Diseq, &di));
331: PetscCall(MatSetValues(H, 1, &zero, 1, &zero, &val, INSERT_VALUES));
332: PetscCall(MatSetValues(H, 1, &one, 1, &one, &two, INSERT_VALUES));
333: }
334: PetscCall(MatAssemblyBegin(H, MAT_FINAL_ASSEMBLY));
335: PetscCall(MatAssemblyEnd(H, MAT_FINAL_ASSEMBLY));
336: if (!user->noeqflag) {
337: PetscCall(VecScatterDestroy(&Descat));
338: PetscCall(VecDestroy(&Deseq));
339: }
340: PetscCall(VecScatterDestroy(&Discat));
341: PetscCall(VecDestroy(&Diseq));
342: PetscFunctionReturn(PETSC_SUCCESS);
343: }
345: /* Evaluate
346: h = [ x0^2 - x1;
347: 1 -(x0^2 - x1)]
348: */
349: PetscErrorCode FormInequalityConstraints(Tao tao, Vec X, Vec CI, void *ctx)
350: {
351: const PetscScalar *x;
352: PetscScalar ci;
353: MPI_Comm comm;
354: PetscMPIInt rank;
355: AppCtx *user = (AppCtx *)ctx;
356: Vec Xseq = user->Xseq;
357: VecScatter scat = user->scat;
359: PetscFunctionBegin;
360: PetscCall(PetscObjectGetComm((PetscObject)tao, &comm));
361: PetscCallMPI(MPI_Comm_rank(comm, &rank));
363: PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD));
364: PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD));
366: if (rank == 0) {
367: PetscCall(VecGetArrayRead(Xseq, &x));
368: ci = x[0] * x[0] - x[1];
369: PetscCall(VecSetValue(CI, 0, ci, INSERT_VALUES));
370: ci = -x[0] * x[0] + x[1] + 1.0;
371: PetscCall(VecSetValue(CI, 1, ci, INSERT_VALUES));
372: PetscCall(VecRestoreArrayRead(Xseq, &x));
373: }
374: PetscCall(VecAssemblyBegin(CI));
375: PetscCall(VecAssemblyEnd(CI));
376: PetscFunctionReturn(PETSC_SUCCESS);
377: }
379: /* Evaluate
380: g = [ x0^2 + x1 - 2]
381: */
382: PetscErrorCode FormEqualityConstraints(Tao tao, Vec X, Vec CE, void *ctx)
383: {
384: const PetscScalar *x;
385: PetscScalar ce;
386: MPI_Comm comm;
387: PetscMPIInt rank;
388: AppCtx *user = (AppCtx *)ctx;
389: Vec Xseq = user->Xseq;
390: VecScatter scat = user->scat;
392: PetscFunctionBegin;
393: PetscCall(PetscObjectGetComm((PetscObject)tao, &comm));
394: PetscCallMPI(MPI_Comm_rank(comm, &rank));
396: PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD));
397: PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD));
399: if (rank == 0) {
400: PetscCall(VecGetArrayRead(Xseq, &x));
401: ce = x[0] * x[0] + x[1] - 2.0;
402: PetscCall(VecSetValue(CE, 0, ce, INSERT_VALUES));
403: PetscCall(VecRestoreArrayRead(Xseq, &x));
404: }
405: PetscCall(VecAssemblyBegin(CE));
406: PetscCall(VecAssemblyEnd(CE));
407: PetscFunctionReturn(PETSC_SUCCESS);
408: }
410: /*
411: grad h = [ 2*x0, -1;
412: -2*x0, 1]
413: */
414: PetscErrorCode FormInequalityJacobian(Tao tao, Vec X, Mat JI, Mat JIpre, void *ctx)
415: {
416: AppCtx *user = (AppCtx *)ctx;
417: PetscInt zero = 0, one = 1, cols[2];
418: PetscScalar vals[2];
419: const PetscScalar *x;
420: Vec Xseq = user->Xseq;
421: VecScatter scat = user->scat;
422: MPI_Comm comm;
423: PetscMPIInt rank;
425: PetscFunctionBegin;
426: PetscCall(PetscObjectGetComm((PetscObject)tao, &comm));
427: PetscCallMPI(MPI_Comm_rank(comm, &rank));
428: PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD));
429: PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD));
431: PetscCall(VecGetArrayRead(Xseq, &x));
432: if (rank == 0) {
433: cols[0] = 0;
434: cols[1] = 1;
435: vals[0] = 2 * x[0];
436: vals[1] = -1.0;
437: PetscCall(MatSetValues(JI, 1, &zero, 2, cols, vals, INSERT_VALUES));
438: vals[0] = -2 * x[0];
439: vals[1] = 1.0;
440: PetscCall(MatSetValues(JI, 1, &one, 2, cols, vals, INSERT_VALUES));
441: }
442: PetscCall(VecRestoreArrayRead(Xseq, &x));
443: PetscCall(MatAssemblyBegin(JI, MAT_FINAL_ASSEMBLY));
444: PetscCall(MatAssemblyEnd(JI, MAT_FINAL_ASSEMBLY));
445: PetscFunctionReturn(PETSC_SUCCESS);
446: }
448: /*
449: grad g = [2*x0, 1.0]
450: */
451: PetscErrorCode FormEqualityJacobian(Tao tao, Vec X, Mat JE, Mat JEpre, void *ctx)
452: {
453: PetscInt zero = 0, cols[2];
454: PetscScalar vals[2];
455: const PetscScalar *x;
456: PetscMPIInt rank;
457: MPI_Comm comm;
459: PetscFunctionBegin;
460: PetscCall(PetscObjectGetComm((PetscObject)tao, &comm));
461: PetscCallMPI(MPI_Comm_rank(comm, &rank));
463: if (rank == 0) {
464: PetscCall(VecGetArrayRead(X, &x));
465: cols[0] = 0;
466: cols[1] = 1;
467: vals[0] = 2 * x[0];
468: vals[1] = 1.0;
469: PetscCall(MatSetValues(JE, 1, &zero, 2, cols, vals, INSERT_VALUES));
470: PetscCall(VecRestoreArrayRead(X, &x));
471: }
472: PetscCall(MatAssemblyBegin(JE, MAT_FINAL_ASSEMBLY));
473: PetscCall(MatAssemblyEnd(JE, MAT_FINAL_ASSEMBLY));
474: PetscFunctionReturn(PETSC_SUCCESS);
475: }
477: /*TEST
479: build:
480: requires: !complex !defined(PETSC_USE_CXX)
482: test:
483: args: -tao_converged_reason -tao_gatol 1.e-6 -tao_type pdipm -tao_pdipm_kkt_shift_pd
484: requires: mumps !single
485: filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"
487: test:
488: suffix: pdipm_2
489: requires: mumps
490: nsize: 2
491: args: -tao_converged_reason -tao_gatol 1.e-6 -tao_type pdipm
492: filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"
494: test:
495: suffix: 2
496: args: -tao_converged_reason
497: filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"
499: test:
500: suffix: 3
501: args: -tao_converged_reason -no_eq
502: filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"
504: test:
505: suffix: 4
506: args: -tao_converged_reason -tao_almm_type classic
507: requires: !single
508: filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"
510: test:
511: suffix: 5
512: args: -tao_converged_reason -tao_almm_type classic -no_eq -tao_almm_subsolver_tao_max_it 100
513: requires: !single !defined(PETSCTEST_VALGRIND)
514: filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"
516: test:
517: suffix: 6
518: args: -tao_converged_reason -tao_almm_subsolver_tao_type bqnktr
519: filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"
521: test:
522: suffix: 7
523: args: -tao_converged_reason -tao_almm_subsolver_tao_type bncg
524: filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"
526: test:
527: suffix: 8
528: nsize: 2
529: args: -tao_converged_reason
530: filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"
532: test:
533: suffix: 9
534: nsize: 2
535: args: -tao_converged_reason -vec_type cuda -mat_type aijcusparse
536: requires: cuda
537: filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"
539: TEST*/