1 /* Program usage: mpiexec -n 2 ./ex1 [-help] [all TAO options] */ 2 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 ---------------------------------------------------------------------- */ 17 18 #include <petsctao.h> 19 20 static char help[] = "Solves constrained optimiztion 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_cmonitor : convergence monitor with constraint norm \n\ 28 -tao_view_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"; 30 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; 44 45 /* -------- User-defined Routines --------- */ 46 PetscErrorCode InitializeProblem(AppCtx *); 47 PetscErrorCode DestroyProblem(AppCtx *); 48 PetscErrorCode FormFunctionGradient(Tao, Vec, PetscReal *, Vec, void *); 49 PetscErrorCode FormHessian(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 *); 54 55 PetscErrorCode main(int argc, char **argv) { 56 Tao tao; 57 KSP ksp; 58 PC pc; 59 AppCtx user; /* application context */ 60 Vec x, G, CI, CE; 61 PetscMPIInt size; 62 TaoType type; 63 PetscReal f; 64 PetscBool pdipm; 65 66 PetscFunctionBeginUser; 67 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 68 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 69 PetscCheck(size <= 2, PETSC_COMM_SELF, PETSC_ERR_WRONG_MPI_SIZE, "More than 2 processors detected. Example written to use max of 2 processors."); 70 71 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "---- Constrained Problem -----\n")); 72 PetscCall(InitializeProblem(&user)); /* sets up problem, function below */ 73 PetscCall(TaoCreate(PETSC_COMM_WORLD, &tao)); 74 PetscCall(TaoSetType(tao, TAOPDIPM)); 75 PetscCall(TaoSetSolution(tao, user.x)); /* gets solution vector from problem */ 76 PetscCall(TaoSetVariableBounds(tao, user.xl, user.xu)); /* sets lower upper bounds from given solution */ 77 PetscCall(TaoSetObjectiveAndGradient(tao, NULL, FormFunctionGradient, (void *)&user)); 78 79 if (!user.noeqflag) PetscCall(TaoSetEqualityConstraintsRoutine(tao, user.ce, FormEqualityConstraints, (void *)&user)); 80 PetscCall(TaoSetInequalityConstraintsRoutine(tao, user.ci, FormInequalityConstraints, (void *)&user)); 81 if (!user.noeqflag) { PetscCall(TaoSetJacobianEqualityRoutine(tao, user.Ae, user.Ae, FormEqualityJacobian, (void *)&user)); /* equality jacobian */ } 82 PetscCall(TaoSetJacobianInequalityRoutine(tao, user.Ai, user.Ai, FormInequalityJacobian, (void *)&user)); /* inequality jacobian */ 83 PetscCall(TaoSetTolerances(tao, 1.e-6, 1.e-6, 1.e-6)); 84 PetscCall(TaoSetConstraintTolerances(tao, 1.e-6, 1.e-6)); 85 86 PetscCall(TaoGetKSP(tao, &ksp)); 87 PetscCall(KSPGetPC(ksp, &pc)); 88 PetscCall(PCSetType(pc, PCCHOLESKY)); 89 /* 90 This algorithm produces matrices with zeros along the diagonal therefore we use 91 MUMPS which provides solver for indefinite matrices 92 */ 93 #if defined(PETSC_HAVE_MUMPS) 94 PetscCall(PCFactorSetMatSolverType(pc, MATSOLVERMUMPS)); /* requires mumps to solve pdipm */ 95 #else 96 PetscCheck(size == 1, PetscObjectComm((PetscObject)pc), PETSC_ERR_SUP, "Requires an external package that supports parallel PCCHOLESKY, e.g., MUMPS."); 97 #endif 98 PetscCall(KSPSetType(ksp, KSPPREONLY)); 99 PetscCall(KSPSetFromOptions(ksp)); 100 101 PetscCall(TaoSetFromOptions(tao)); 102 PetscCall(TaoGetType(tao, &type)); 103 PetscCall(PetscObjectTypeCompare((PetscObject)tao, TAOPDIPM, &pdipm)); 104 if (pdipm) { 105 PetscCall(TaoSetHessian(tao, user.H, user.H, FormHessian, (void *)&user)); 106 if (user.initview) { 107 PetscCall(TaoSetUp(tao)); 108 PetscCall(VecDuplicate(user.x, &G)); 109 PetscCall(FormFunctionGradient(tao, user.x, &f, G, (void *)&user)); 110 PetscCall(FormHessian(tao, user.x, user.H, user.H, (void *)&user)); 111 PetscCall(PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_WORLD)); 112 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial point X:\n")); 113 PetscCall(VecView(user.x, PETSC_VIEWER_STDOUT_WORLD)); 114 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial objective f(x) = %g\n", (double)f)); 115 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial gradient and Hessian:\n")); 116 PetscCall(VecView(G, PETSC_VIEWER_STDOUT_WORLD)); 117 PetscCall(MatView(user.H, PETSC_VIEWER_STDOUT_WORLD)); 118 PetscCall(VecDestroy(&G)); 119 PetscCall(FormInequalityJacobian(tao, user.x, user.Ai, user.Ai, (void *)&user)); 120 PetscCall(MatCreateVecs(user.Ai, NULL, &CI)); 121 PetscCall(FormInequalityConstraints(tao, user.x, CI, (void *)&user)); 122 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial inequality constraints and Jacobian:\n")); 123 PetscCall(VecView(CI, PETSC_VIEWER_STDOUT_WORLD)); 124 PetscCall(MatView(user.Ai, PETSC_VIEWER_STDOUT_WORLD)); 125 PetscCall(VecDestroy(&CI)); 126 if (!user.noeqflag) { 127 PetscCall(FormEqualityJacobian(tao, user.x, user.Ae, user.Ae, (void *)&user)); 128 PetscCall(MatCreateVecs(user.Ae, NULL, &CE)); 129 PetscCall(FormEqualityConstraints(tao, user.x, CE, (void *)&user)); 130 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial equality constraints and Jacobian:\n")); 131 PetscCall(VecView(CE, PETSC_VIEWER_STDOUT_WORLD)); 132 PetscCall(MatView(user.Ae, PETSC_VIEWER_STDOUT_WORLD)); 133 PetscCall(VecDestroy(&CE)); 134 } 135 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\n")); 136 PetscCall(PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_WORLD)); 137 } 138 } 139 140 PetscCall(TaoSolve(tao)); 141 PetscCall(TaoGetSolution(tao, &x)); 142 PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD)); 143 144 /* Free objects */ 145 PetscCall(DestroyProblem(&user)); 146 PetscCall(TaoDestroy(&tao)); 147 PetscCall(PetscFinalize()); 148 return 0; 149 } 150 151 PetscErrorCode InitializeProblem(AppCtx *user) { 152 PetscMPIInt size; 153 PetscMPIInt rank; 154 PetscInt nloc, neloc, niloc; 155 156 PetscFunctionBegin; 157 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 158 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 159 user->noeqflag = PETSC_FALSE; 160 PetscCall(PetscOptionsGetBool(NULL, NULL, "-no_eq", &user->noeqflag, NULL)); 161 user->initview = PETSC_FALSE; 162 PetscCall(PetscOptionsGetBool(NULL, NULL, "-init_view", &user->initview, NULL)); 163 164 if (!user->noeqflag) { 165 /* Tell user the correct solution, not an error checking */ 166 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Solution should be f(1,1)=-2\n")); 167 } 168 169 /* create vector x and set initial values */ 170 user->n = 2; /* global length */ 171 nloc = (size == 1) ? user->n : 1; 172 PetscCall(VecCreate(PETSC_COMM_WORLD, &user->x)); 173 PetscCall(VecSetSizes(user->x, nloc, user->n)); 174 PetscCall(VecSetFromOptions(user->x)); 175 PetscCall(VecSet(user->x, 0)); 176 177 /* create and set lower and upper bound vectors */ 178 PetscCall(VecDuplicate(user->x, &user->xl)); 179 PetscCall(VecDuplicate(user->x, &user->xu)); 180 PetscCall(VecSet(user->xl, -1.0)); 181 PetscCall(VecSet(user->xu, 2.0)); 182 183 /* create scater to zero */ 184 PetscCall(VecScatterCreateToZero(user->x, &user->scat, &user->Xseq)); 185 186 user->ne = 1; 187 user->ni = 2; 188 neloc = (rank == 0) ? user->ne : 0; 189 niloc = (size == 1) ? user->ni : 1; 190 191 if (!user->noeqflag) { 192 PetscCall(VecCreate(PETSC_COMM_WORLD, &user->ce)); /* a 1x1 vec for equality constraints */ 193 PetscCall(VecSetSizes(user->ce, neloc, user->ne)); 194 PetscCall(VecSetFromOptions(user->ce)); 195 PetscCall(VecSetUp(user->ce)); 196 } 197 198 PetscCall(VecCreate(PETSC_COMM_WORLD, &user->ci)); /* a 2x1 vec for inequality constraints */ 199 PetscCall(VecSetSizes(user->ci, niloc, user->ni)); 200 PetscCall(VecSetFromOptions(user->ci)); 201 PetscCall(VecSetUp(user->ci)); 202 203 /* nexn & nixn matricies for equally and inequalty constraints */ 204 if (!user->noeqflag) { 205 PetscCall(MatCreate(PETSC_COMM_WORLD, &user->Ae)); 206 PetscCall(MatSetSizes(user->Ae, neloc, nloc, user->ne, user->n)); 207 PetscCall(MatSetFromOptions(user->Ae)); 208 PetscCall(MatSetUp(user->Ae)); 209 } 210 211 PetscCall(MatCreate(PETSC_COMM_WORLD, &user->Ai)); 212 PetscCall(MatSetSizes(user->Ai, niloc, nloc, user->ni, user->n)); 213 PetscCall(MatSetFromOptions(user->Ai)); 214 PetscCall(MatSetUp(user->Ai)); 215 216 PetscCall(MatCreate(PETSC_COMM_WORLD, &user->H)); 217 PetscCall(MatSetSizes(user->H, nloc, nloc, user->n, user->n)); 218 PetscCall(MatSetFromOptions(user->H)); 219 PetscCall(MatSetUp(user->H)); 220 PetscFunctionReturn(0); 221 } 222 223 PetscErrorCode DestroyProblem(AppCtx *user) { 224 PetscFunctionBegin; 225 if (!user->noeqflag) PetscCall(MatDestroy(&user->Ae)); 226 PetscCall(MatDestroy(&user->Ai)); 227 PetscCall(MatDestroy(&user->H)); 228 229 PetscCall(VecDestroy(&user->x)); 230 if (!user->noeqflag) PetscCall(VecDestroy(&user->ce)); 231 PetscCall(VecDestroy(&user->ci)); 232 PetscCall(VecDestroy(&user->xl)); 233 PetscCall(VecDestroy(&user->xu)); 234 PetscCall(VecDestroy(&user->Xseq)); 235 PetscCall(VecScatterDestroy(&user->scat)); 236 PetscFunctionReturn(0); 237 } 238 239 /* Evaluate 240 f(x) = (x0 - 2)^2 + (x1 - 2)^2 - 2*(x0 + x1) 241 G = grad f = [2*(x0 - 2) - 2; 242 2*(x1 - 2) - 2] 243 */ 244 PetscErrorCode FormFunctionGradient(Tao tao, Vec X, PetscReal *f, Vec G, void *ctx) { 245 PetscScalar g; 246 const PetscScalar *x; 247 MPI_Comm comm; 248 PetscMPIInt rank; 249 PetscReal fin; 250 AppCtx *user = (AppCtx *)ctx; 251 Vec Xseq = user->Xseq; 252 VecScatter scat = user->scat; 253 254 PetscFunctionBegin; 255 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 256 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 257 258 PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 259 PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 260 261 fin = 0.0; 262 if (rank == 0) { 263 PetscCall(VecGetArrayRead(Xseq, &x)); 264 fin = (x[0] - 2.0) * (x[0] - 2.0) + (x[1] - 2.0) * (x[1] - 2.0) - 2.0 * (x[0] + x[1]); 265 g = 2.0 * (x[0] - 2.0) - 2.0; 266 PetscCall(VecSetValue(G, 0, g, INSERT_VALUES)); 267 g = 2.0 * (x[1] - 2.0) - 2.0; 268 PetscCall(VecSetValue(G, 1, g, INSERT_VALUES)); 269 PetscCall(VecRestoreArrayRead(Xseq, &x)); 270 } 271 PetscCallMPI(MPI_Allreduce(&fin, f, 1, MPIU_REAL, MPIU_SUM, comm)); 272 PetscCall(VecAssemblyBegin(G)); 273 PetscCall(VecAssemblyEnd(G)); 274 PetscFunctionReturn(0); 275 } 276 277 /* Evaluate 278 H = fxx + grad (grad g^T*DI) - grad (grad h^T*DE)] 279 = [ 2*(1+de[0]-di[0]+di[1]), 0; 280 0, 2] 281 */ 282 PetscErrorCode FormHessian(Tao tao, Vec x, Mat H, Mat Hpre, void *ctx) { 283 AppCtx *user = (AppCtx *)ctx; 284 Vec DE, DI; 285 const PetscScalar *de, *di; 286 PetscInt zero = 0, one = 1; 287 PetscScalar two = 2.0; 288 PetscScalar val = 0.0; 289 Vec Deseq, Diseq; 290 VecScatter Descat, Discat; 291 PetscMPIInt rank; 292 MPI_Comm comm; 293 294 PetscFunctionBegin; 295 PetscCall(TaoGetDualVariables(tao, &DE, &DI)); 296 297 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 298 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 299 300 if (!user->noeqflag) { 301 PetscCall(VecScatterCreateToZero(DE, &Descat, &Deseq)); 302 PetscCall(VecScatterBegin(Descat, DE, Deseq, INSERT_VALUES, SCATTER_FORWARD)); 303 PetscCall(VecScatterEnd(Descat, DE, Deseq, INSERT_VALUES, SCATTER_FORWARD)); 304 } 305 PetscCall(VecScatterCreateToZero(DI, &Discat, &Diseq)); 306 PetscCall(VecScatterBegin(Discat, DI, Diseq, INSERT_VALUES, SCATTER_FORWARD)); 307 PetscCall(VecScatterEnd(Discat, DI, Diseq, INSERT_VALUES, SCATTER_FORWARD)); 308 309 if (rank == 0) { 310 if (!user->noeqflag) { PetscCall(VecGetArrayRead(Deseq, &de)); /* places equality constraint dual into array */ } 311 PetscCall(VecGetArrayRead(Diseq, &di)); /* places inequality constraint dual into array */ 312 313 if (!user->noeqflag) { 314 val = 2.0 * (1 + de[0] - di[0] + di[1]); 315 PetscCall(VecRestoreArrayRead(Deseq, &de)); 316 PetscCall(VecRestoreArrayRead(Diseq, &di)); 317 } else { 318 val = 2.0 * (1 - di[0] + di[1]); 319 } 320 PetscCall(VecRestoreArrayRead(Diseq, &di)); 321 PetscCall(MatSetValues(H, 1, &zero, 1, &zero, &val, INSERT_VALUES)); 322 PetscCall(MatSetValues(H, 1, &one, 1, &one, &two, INSERT_VALUES)); 323 } 324 PetscCall(MatAssemblyBegin(H, MAT_FINAL_ASSEMBLY)); 325 PetscCall(MatAssemblyEnd(H, MAT_FINAL_ASSEMBLY)); 326 if (!user->noeqflag) { 327 PetscCall(VecScatterDestroy(&Descat)); 328 PetscCall(VecDestroy(&Deseq)); 329 } 330 PetscCall(VecScatterDestroy(&Discat)); 331 PetscCall(VecDestroy(&Diseq)); 332 PetscFunctionReturn(0); 333 } 334 335 /* Evaluate 336 h = [ x0^2 - x1; 337 1 -(x0^2 - x1)] 338 */ 339 PetscErrorCode FormInequalityConstraints(Tao tao, Vec X, Vec CI, void *ctx) { 340 const PetscScalar *x; 341 PetscScalar ci; 342 MPI_Comm comm; 343 PetscMPIInt rank; 344 AppCtx *user = (AppCtx *)ctx; 345 Vec Xseq = user->Xseq; 346 VecScatter scat = user->scat; 347 348 PetscFunctionBegin; 349 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 350 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 351 352 PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 353 PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 354 355 if (rank == 0) { 356 PetscCall(VecGetArrayRead(Xseq, &x)); 357 ci = x[0] * x[0] - x[1]; 358 PetscCall(VecSetValue(CI, 0, ci, INSERT_VALUES)); 359 ci = -x[0] * x[0] + x[1] + 1.0; 360 PetscCall(VecSetValue(CI, 1, ci, INSERT_VALUES)); 361 PetscCall(VecRestoreArrayRead(Xseq, &x)); 362 } 363 PetscCall(VecAssemblyBegin(CI)); 364 PetscCall(VecAssemblyEnd(CI)); 365 PetscFunctionReturn(0); 366 } 367 368 /* Evaluate 369 g = [ x0^2 + x1 - 2] 370 */ 371 PetscErrorCode FormEqualityConstraints(Tao tao, Vec X, Vec CE, void *ctx) { 372 const PetscScalar *x; 373 PetscScalar ce; 374 MPI_Comm comm; 375 PetscMPIInt rank; 376 AppCtx *user = (AppCtx *)ctx; 377 Vec Xseq = user->Xseq; 378 VecScatter scat = user->scat; 379 380 PetscFunctionBegin; 381 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 382 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 383 384 PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 385 PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 386 387 if (rank == 0) { 388 PetscCall(VecGetArrayRead(Xseq, &x)); 389 ce = x[0] * x[0] + x[1] - 2.0; 390 PetscCall(VecSetValue(CE, 0, ce, INSERT_VALUES)); 391 PetscCall(VecRestoreArrayRead(Xseq, &x)); 392 } 393 PetscCall(VecAssemblyBegin(CE)); 394 PetscCall(VecAssemblyEnd(CE)); 395 PetscFunctionReturn(0); 396 } 397 398 /* 399 grad h = [ 2*x0, -1; 400 -2*x0, 1] 401 */ 402 PetscErrorCode FormInequalityJacobian(Tao tao, Vec X, Mat JI, Mat JIpre, void *ctx) { 403 AppCtx *user = (AppCtx *)ctx; 404 PetscInt zero = 0, one = 1, cols[2]; 405 PetscScalar vals[2]; 406 const PetscScalar *x; 407 Vec Xseq = user->Xseq; 408 VecScatter scat = user->scat; 409 MPI_Comm comm; 410 PetscMPIInt rank; 411 412 PetscFunctionBegin; 413 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 414 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 415 PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 416 PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 417 418 PetscCall(VecGetArrayRead(Xseq, &x)); 419 if (rank == 0) { 420 cols[0] = 0; 421 cols[1] = 1; 422 vals[0] = 2 * x[0]; 423 vals[1] = -1.0; 424 PetscCall(MatSetValues(JI, 1, &zero, 2, cols, vals, INSERT_VALUES)); 425 vals[0] = -2 * x[0]; 426 vals[1] = 1.0; 427 PetscCall(MatSetValues(JI, 1, &one, 2, cols, vals, INSERT_VALUES)); 428 } 429 PetscCall(VecRestoreArrayRead(Xseq, &x)); 430 PetscCall(MatAssemblyBegin(JI, MAT_FINAL_ASSEMBLY)); 431 PetscCall(MatAssemblyEnd(JI, MAT_FINAL_ASSEMBLY)); 432 PetscFunctionReturn(0); 433 } 434 435 /* 436 grad g = [2*x0 437 1.0 ] 438 */ 439 PetscErrorCode FormEqualityJacobian(Tao tao, Vec X, Mat JE, Mat JEpre, void *ctx) { 440 PetscInt zero = 0, cols[2]; 441 PetscScalar vals[2]; 442 const PetscScalar *x; 443 PetscMPIInt rank; 444 MPI_Comm comm; 445 446 PetscFunctionBegin; 447 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 448 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 449 450 if (rank == 0) { 451 PetscCall(VecGetArrayRead(X, &x)); 452 cols[0] = 0; 453 cols[1] = 1; 454 vals[0] = 2 * x[0]; 455 vals[1] = 1.0; 456 PetscCall(MatSetValues(JE, 1, &zero, 2, cols, vals, INSERT_VALUES)); 457 PetscCall(VecRestoreArrayRead(X, &x)); 458 } 459 PetscCall(MatAssemblyBegin(JE, MAT_FINAL_ASSEMBLY)); 460 PetscCall(MatAssemblyEnd(JE, MAT_FINAL_ASSEMBLY)); 461 PetscFunctionReturn(0); 462 } 463 464 /*TEST 465 466 build: 467 requires: !complex !defined(PETSC_USE_CXX) 468 469 test: 470 args: -tao_converged_reason -tao_pdipm_kkt_shift_pd 471 requires: mumps 472 473 test: 474 suffix: 2 475 nsize: 2 476 args: -tao_converged_reason -tao_pdipm_kkt_shift_pd 477 requires: mumps 478 479 test: 480 suffix: 3 481 args: -tao_converged_reason -no_eq 482 requires: mumps 483 484 test: 485 suffix: 4 486 nsize: 2 487 args: -tao_converged_reason -no_eq 488 requires: mumps 489 490 test: 491 suffix: 5 492 args: -tao_cmonitor -tao_type almm 493 requires: mumps 494 495 test: 496 suffix: 6 497 args: -tao_cmonitor -tao_type almm -tao_almm_type phr 498 requires: mumps 499 500 test: 501 suffix: 7 502 nsize: 2 503 requires: mumps 504 args: -tao_cmonitor -tao_type almm 505 506 test: 507 suffix: 8 508 nsize: 2 509 requires: cuda mumps 510 args: -tao_cmonitor -tao_type almm -vec_type cuda -mat_type aijcusparse 511 512 test: 513 suffix: 9 514 nsize: 2 515 args: -tao_cmonitor -tao_type almm -no_eq 516 requires: mumps 517 518 test: 519 suffix: 10 520 nsize: 2 521 args: -tao_cmonitor -tao_type almm -tao_almm_type phr -no_eq 522 requires: mumps 523 524 TEST*/ 525