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