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