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, mumps; 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_WORLD, 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 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)); 84 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)); 91 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(KSPSetFromOptions(ksp)); 106 PetscCall(TaoSetHessian(tao, user.H, user.H, FormHessian, (void *)&user)); 107 } 108 109 /* Print out an initial view of the problem */ 110 if (user.initview) { 111 PetscCall(TaoSetUp(tao)); 112 PetscCall(VecDuplicate(user.x, &G)); 113 PetscCall(FormFunctionGradient(tao, user.x, &f, G, (void *)&user)); 114 PetscCall(FormHessian(tao, user.x, user.H, user.H, (void *)&user)); 115 PetscCall(PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_WORLD)); 116 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial point X:\n")); 117 PetscCall(VecView(user.x, PETSC_VIEWER_STDOUT_WORLD)); 118 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial objective f(x) = %g\n", (double)f)); 119 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial gradient and Hessian:\n")); 120 PetscCall(VecView(G, PETSC_VIEWER_STDOUT_WORLD)); 121 PetscCall(MatView(user.H, PETSC_VIEWER_STDOUT_WORLD)); 122 PetscCall(VecDestroy(&G)); 123 PetscCall(FormInequalityJacobian(tao, user.x, user.Ai, user.Ai, (void *)&user)); 124 PetscCall(MatCreateVecs(user.Ai, NULL, &CI)); 125 PetscCall(FormInequalityConstraints(tao, user.x, CI, (void *)&user)); 126 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial inequality constraints and Jacobian:\n")); 127 PetscCall(VecView(CI, PETSC_VIEWER_STDOUT_WORLD)); 128 PetscCall(MatView(user.Ai, PETSC_VIEWER_STDOUT_WORLD)); 129 PetscCall(VecDestroy(&CI)); 130 if (!user.noeqflag) { 131 PetscCall(FormEqualityJacobian(tao, user.x, user.Ae, user.Ae, (void *)&user)); 132 PetscCall(MatCreateVecs(user.Ae, NULL, &CE)); 133 PetscCall(FormEqualityConstraints(tao, user.x, CE, (void *)&user)); 134 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\nInitial equality constraints and Jacobian:\n")); 135 PetscCall(VecView(CE, PETSC_VIEWER_STDOUT_WORLD)); 136 PetscCall(MatView(user.Ae, PETSC_VIEWER_STDOUT_WORLD)); 137 PetscCall(VecDestroy(&CE)); 138 } 139 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\n")); 140 PetscCall(PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_WORLD)); 141 } 142 143 PetscCall(TaoSolve(tao)); 144 PetscCall(TaoGetSolution(tao, &x)); 145 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Found solution:\n")); 146 PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD)); 147 148 /* Free objects */ 149 PetscCall(DestroyProblem(&user)); 150 PetscCall(TaoDestroy(&tao)); 151 PetscCall(PetscFinalize()); 152 return 0; 153 } 154 155 PetscErrorCode InitializeProblem(AppCtx *user) 156 { 157 PetscMPIInt size; 158 PetscMPIInt rank; 159 PetscInt nloc, neloc, niloc; 160 161 PetscFunctionBegin; 162 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 163 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 164 user->noeqflag = PETSC_FALSE; 165 PetscCall(PetscOptionsGetBool(NULL, NULL, "-no_eq", &user->noeqflag, NULL)); 166 user->initview = PETSC_FALSE; 167 PetscCall(PetscOptionsGetBool(NULL, NULL, "-init_view", &user->initview, NULL)); 168 169 /* Tell user the correct solution, not an error checking */ 170 if (!user->noeqflag) { 171 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Expected solution: f(1, 1) = -2\n")); 172 } else { 173 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Expected solution (-no_eq): f(1.73205, 2) = -7.3923\n")); 174 } 175 176 /* create vector x and set initial values */ 177 user->n = 2; /* global length */ 178 nloc = (size == 1) ? user->n : 1; 179 PetscCall(VecCreate(PETSC_COMM_WORLD, &user->x)); 180 PetscCall(VecSetSizes(user->x, nloc, user->n)); 181 PetscCall(VecSetFromOptions(user->x)); 182 PetscCall(VecSet(user->x, 0.0)); 183 184 /* create and set lower and upper bound vectors */ 185 PetscCall(VecDuplicate(user->x, &user->xl)); 186 PetscCall(VecDuplicate(user->x, &user->xu)); 187 PetscCall(VecSet(user->xl, -1.0)); 188 PetscCall(VecSet(user->xu, 2.0)); 189 190 /* create scater to zero */ 191 PetscCall(VecScatterCreateToZero(user->x, &user->scat, &user->Xseq)); 192 193 user->ne = 1; 194 user->ni = 2; 195 neloc = (rank == 0) ? user->ne : 0; 196 niloc = (size == 1) ? user->ni : 1; 197 198 if (!user->noeqflag) { 199 PetscCall(VecCreate(PETSC_COMM_WORLD, &user->ce)); /* a 1x1 vec for equality constraints */ 200 PetscCall(VecSetSizes(user->ce, neloc, user->ne)); 201 PetscCall(VecSetFromOptions(user->ce)); 202 PetscCall(VecSetUp(user->ce)); 203 } 204 205 PetscCall(VecCreate(PETSC_COMM_WORLD, &user->ci)); /* a 2x1 vec for inequality constraints */ 206 PetscCall(VecSetSizes(user->ci, niloc, user->ni)); 207 PetscCall(VecSetFromOptions(user->ci)); 208 PetscCall(VecSetUp(user->ci)); 209 210 /* nexn & nixn matricies for equally and inequalty constraints */ 211 if (!user->noeqflag) { 212 PetscCall(MatCreate(PETSC_COMM_WORLD, &user->Ae)); 213 PetscCall(MatSetSizes(user->Ae, neloc, nloc, user->ne, user->n)); 214 PetscCall(MatSetFromOptions(user->Ae)); 215 PetscCall(MatSetUp(user->Ae)); 216 } 217 218 PetscCall(MatCreate(PETSC_COMM_WORLD, &user->Ai)); 219 PetscCall(MatSetSizes(user->Ai, niloc, nloc, user->ni, user->n)); 220 PetscCall(MatSetFromOptions(user->Ai)); 221 PetscCall(MatSetUp(user->Ai)); 222 223 PetscCall(MatCreate(PETSC_COMM_WORLD, &user->H)); 224 PetscCall(MatSetSizes(user->H, nloc, nloc, user->n, user->n)); 225 PetscCall(MatSetFromOptions(user->H)); 226 PetscCall(MatSetUp(user->H)); 227 PetscFunctionReturn(0); 228 } 229 230 PetscErrorCode DestroyProblem(AppCtx *user) 231 { 232 PetscFunctionBegin; 233 if (!user->noeqflag) PetscCall(MatDestroy(&user->Ae)); 234 PetscCall(MatDestroy(&user->Ai)); 235 PetscCall(MatDestroy(&user->H)); 236 237 PetscCall(VecDestroy(&user->x)); 238 if (!user->noeqflag) PetscCall(VecDestroy(&user->ce)); 239 PetscCall(VecDestroy(&user->ci)); 240 PetscCall(VecDestroy(&user->xl)); 241 PetscCall(VecDestroy(&user->xu)); 242 PetscCall(VecDestroy(&user->Xseq)); 243 PetscCall(VecScatterDestroy(&user->scat)); 244 PetscFunctionReturn(0); 245 } 246 247 /* Evaluate 248 f(x) = (x0 - 2)^2 + (x1 - 2)^2 - 2*(x0 + x1) 249 G = grad f = [2*(x0 - 2) - 2; 250 2*(x1 - 2) - 2] 251 */ 252 PetscErrorCode FormFunctionGradient(Tao tao, Vec X, PetscReal *f, Vec G, void *ctx) 253 { 254 PetscScalar g; 255 const PetscScalar *x; 256 MPI_Comm comm; 257 PetscMPIInt rank; 258 PetscReal fin; 259 AppCtx *user = (AppCtx *)ctx; 260 Vec Xseq = user->Xseq; 261 VecScatter scat = user->scat; 262 263 PetscFunctionBegin; 264 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 265 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 266 267 PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 268 PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 269 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 g = 2.0 * (x[0] - 2.0) - 2.0; 275 PetscCall(VecSetValue(G, 0, g, INSERT_VALUES)); 276 g = 2.0 * (x[1] - 2.0) - 2.0; 277 PetscCall(VecSetValue(G, 1, g, INSERT_VALUES)); 278 PetscCall(VecRestoreArrayRead(Xseq, &x)); 279 } 280 PetscCallMPI(MPI_Allreduce(&fin, f, 1, MPIU_REAL, MPIU_SUM, comm)); 281 PetscCall(VecAssemblyBegin(G)); 282 PetscCall(VecAssemblyEnd(G)); 283 PetscFunctionReturn(0); 284 } 285 286 /* Evaluate 287 H = fxx + grad (grad g^T*DI) - grad (grad h^T*DE)] 288 = [ 2*(1+de[0]-di[0]+di[1]), 0; 289 0, 2] 290 */ 291 PetscErrorCode FormHessian(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; 303 304 PetscFunctionBegin; 305 PetscCall(TaoGetDualVariables(tao, &DE, &DI)); 306 307 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 308 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 309 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)); 318 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 */ 322 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(0); 343 } 344 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; 358 359 PetscFunctionBegin; 360 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 361 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 362 363 PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 364 PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 365 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(0); 377 } 378 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; 391 392 PetscFunctionBegin; 393 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 394 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 395 396 PetscCall(VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 397 PetscCall(VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD)); 398 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(0); 408 } 409 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; 424 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)); 430 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(0); 446 } 447 448 /* 449 grad g = [2*x0 450 1.0 ] 451 */ 452 PetscErrorCode FormEqualityJacobian(Tao tao, Vec X, Mat JE, Mat JEpre, void *ctx) 453 { 454 PetscInt zero = 0, cols[2]; 455 PetscScalar vals[2]; 456 const PetscScalar *x; 457 PetscMPIInt rank; 458 MPI_Comm comm; 459 460 PetscFunctionBegin; 461 PetscCall(PetscObjectGetComm((PetscObject)tao, &comm)); 462 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 463 464 if (rank == 0) { 465 PetscCall(VecGetArrayRead(X, &x)); 466 cols[0] = 0; 467 cols[1] = 1; 468 vals[0] = 2 * x[0]; 469 vals[1] = 1.0; 470 PetscCall(MatSetValues(JE, 1, &zero, 2, cols, vals, INSERT_VALUES)); 471 PetscCall(VecRestoreArrayRead(X, &x)); 472 } 473 PetscCall(MatAssemblyBegin(JE, MAT_FINAL_ASSEMBLY)); 474 PetscCall(MatAssemblyEnd(JE, MAT_FINAL_ASSEMBLY)); 475 PetscFunctionReturn(0); 476 } 477 478 /*TEST 479 480 build: 481 requires: !complex !defined(PETSC_USE_CXX) 482 483 test: 484 args: -tao_converged_reason -tao_gatol 1.e-6 -tao_type pdipm -tao_pdipm_kkt_shift_pd 485 requires: mumps 486 filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g" 487 488 test: 489 suffix: 2 490 args: -tao_converged_reason 491 filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g" 492 493 test: 494 suffix: 3 495 args: -tao_converged_reason -no_eq 496 filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g" 497 498 test: 499 suffix: 4 500 args: -tao_converged_reason -tao_almm_type classic 501 requires: !single 502 filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g" 503 504 test: 505 suffix: 5 506 args: -tao_converged_reason -tao_almm_type classic -no_eq 507 requires: !single 508 filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g" 509 510 test: 511 suffix: 6 512 args: -tao_converged_reason -tao_almm_subsolver_tao_type bqnktr 513 filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g" 514 515 test: 516 suffix: 7 517 args: -tao_converged_reason -tao_almm_subsolver_tao_type bncg 518 filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g" 519 520 test: 521 suffix: 8 522 nsize: 2 523 args: -tao_converged_reason 524 filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g" 525 526 test: 527 suffix: 9 528 nsize: 2 529 args: -tao_converged_reason -vec_type cuda -mat_type aijcusparse 530 requires: cuda 531 filter: sed -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g" 532 533 TEST*/ 534