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 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";

 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 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 *);

 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;

 68:   PetscInitialize(&argc, &argv, (char *)0, help);
 69:   MPI_Comm_size(PETSC_COMM_WORLD, &size);

 72:   PetscPrintf(PETSC_COMM_WORLD, "---- Constrained Problem -----\n");
 73:   InitializeProblem(&user); /* sets up problem, function below */

 75:   TaoCreate(PETSC_COMM_WORLD, &tao);
 76:   TaoSetType(tao, TAOALMM);
 77:   TaoSetSolution(tao, user.x);
 78:   TaoSetVariableBounds(tao, user.xl, user.xu);
 79:   TaoSetObjectiveAndGradient(tao, NULL, FormFunctionGradient, (void *)&user);
 80:   TaoSetTolerances(tao, 1.e-4, 0.0, 0.0);
 81:   TaoSetConstraintTolerances(tao, 1.e-4, 0.0);
 82:   TaoSetMaximumFunctionEvaluations(tao, 1e6);
 83:   TaoSetFromOptions(tao);

 85:   if (!user.noeqflag) {
 86:     TaoSetEqualityConstraintsRoutine(tao, user.ce, FormEqualityConstraints, (void *)&user);
 87:     TaoSetJacobianEqualityRoutine(tao, user.Ae, user.Ae, FormEqualityJacobian, (void *)&user);
 88:   }
 89:   TaoSetInequalityConstraintsRoutine(tao, user.ci, FormInequalityConstraints, (void *)&user);
 90:   TaoSetJacobianInequalityRoutine(tao, user.Ai, user.Ai, FormInequalityJacobian, (void *)&user);

 92:   TaoGetType(tao, &type);
 93:   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:     PetscHasExternalPackage("mumps", &mumps);
101:     TaoGetKSP(tao, &ksp);
102:     KSPGetPC(ksp, &pc);
103:     PCFactorSetMatSolverType(pc, MATSOLVERMUMPS);
104:     KSPSetType(ksp, KSPPREONLY);
105:     KSPSetFromOptions(ksp);
106:     TaoSetHessian(tao, user.H, user.H, FormHessian, (void *)&user);
107:   }

109:   /* Print out an initial view of the problem */
110:   if (user.initview) {
111:     TaoSetUp(tao);
112:     VecDuplicate(user.x, &G);
113:     FormFunctionGradient(tao, user.x, &f, G, (void *)&user);
114:     FormHessian(tao, user.x, user.H, user.H, (void *)&user);
115:     PetscViewerASCIIPushTab(PETSC_VIEWER_STDOUT_WORLD);
116:     PetscPrintf(PETSC_COMM_WORLD, "\nInitial point X:\n");
117:     VecView(user.x, PETSC_VIEWER_STDOUT_WORLD);
118:     PetscPrintf(PETSC_COMM_WORLD, "\nInitial objective f(x) = %g\n", (double)f);
119:     PetscPrintf(PETSC_COMM_WORLD, "\nInitial gradient and Hessian:\n");
120:     VecView(G, PETSC_VIEWER_STDOUT_WORLD);
121:     MatView(user.H, PETSC_VIEWER_STDOUT_WORLD);
122:     VecDestroy(&G);
123:     FormInequalityJacobian(tao, user.x, user.Ai, user.Ai, (void *)&user);
124:     MatCreateVecs(user.Ai, NULL, &CI);
125:     FormInequalityConstraints(tao, user.x, CI, (void *)&user);
126:     PetscPrintf(PETSC_COMM_WORLD, "\nInitial inequality constraints and Jacobian:\n");
127:     VecView(CI, PETSC_VIEWER_STDOUT_WORLD);
128:     MatView(user.Ai, PETSC_VIEWER_STDOUT_WORLD);
129:     VecDestroy(&CI);
130:     if (!user.noeqflag) {
131:       FormEqualityJacobian(tao, user.x, user.Ae, user.Ae, (void *)&user);
132:       MatCreateVecs(user.Ae, NULL, &CE);
133:       FormEqualityConstraints(tao, user.x, CE, (void *)&user);
134:       PetscPrintf(PETSC_COMM_WORLD, "\nInitial equality constraints and Jacobian:\n");
135:       VecView(CE, PETSC_VIEWER_STDOUT_WORLD);
136:       MatView(user.Ae, PETSC_VIEWER_STDOUT_WORLD);
137:       VecDestroy(&CE);
138:     }
139:     PetscPrintf(PETSC_COMM_WORLD, "\n");
140:     PetscViewerASCIIPopTab(PETSC_VIEWER_STDOUT_WORLD);
141:   }

143:   TaoSolve(tao);
144:   TaoGetSolution(tao, &x);
145:   PetscPrintf(PETSC_COMM_WORLD, "Found solution:\n");
146:   VecView(x, PETSC_VIEWER_STDOUT_WORLD);

148:   /* Free objects */
149:   DestroyProblem(&user);
150:   TaoDestroy(&tao);
151:   PetscFinalize();
152:   return 0;
153: }

155: PetscErrorCode InitializeProblem(AppCtx *user)
156: {
157:   PetscMPIInt size;
158:   PetscMPIInt rank;
159:   PetscInt    nloc, neloc, niloc;

161:   MPI_Comm_size(PETSC_COMM_WORLD, &size);
162:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
163:   user->noeqflag = PETSC_FALSE;
164:   PetscOptionsGetBool(NULL, NULL, "-no_eq", &user->noeqflag, NULL);
165:   user->initview = PETSC_FALSE;
166:   PetscOptionsGetBool(NULL, NULL, "-init_view", &user->initview, NULL);

168:   /* Tell user the correct solution, not an error checking */
169:   if (!user->noeqflag) {
170:     PetscPrintf(PETSC_COMM_WORLD, "Expected solution: f(1, 1) = -2\n");
171:   } else {
172:     PetscPrintf(PETSC_COMM_WORLD, "Expected solution (-no_eq): f(1.73205, 2) = -7.3923\n");
173:   }

175:   /* create vector x and set initial values */
176:   user->n = 2; /* global length */
177:   nloc    = (size == 1) ? user->n : 1;
178:   VecCreate(PETSC_COMM_WORLD, &user->x);
179:   VecSetSizes(user->x, nloc, user->n);
180:   VecSetFromOptions(user->x);
181:   VecSet(user->x, 0.0);

183:   /* create and set lower and upper bound vectors */
184:   VecDuplicate(user->x, &user->xl);
185:   VecDuplicate(user->x, &user->xu);
186:   VecSet(user->xl, -1.0);
187:   VecSet(user->xu, 2.0);

189:   /* create scater to zero */
190:   VecScatterCreateToZero(user->x, &user->scat, &user->Xseq);

192:   user->ne = 1;
193:   user->ni = 2;
194:   neloc    = (rank == 0) ? user->ne : 0;
195:   niloc    = (size == 1) ? user->ni : 1;

197:   if (!user->noeqflag) {
198:     VecCreate(PETSC_COMM_WORLD, &user->ce); /* a 1x1 vec for equality constraints */
199:     VecSetSizes(user->ce, neloc, user->ne);
200:     VecSetFromOptions(user->ce);
201:     VecSetUp(user->ce);
202:   }

204:   VecCreate(PETSC_COMM_WORLD, &user->ci); /* a 2x1 vec for inequality constraints */
205:   VecSetSizes(user->ci, niloc, user->ni);
206:   VecSetFromOptions(user->ci);
207:   VecSetUp(user->ci);

209:   /* nexn & nixn matricies for equally and inequalty constraints */
210:   if (!user->noeqflag) {
211:     MatCreate(PETSC_COMM_WORLD, &user->Ae);
212:     MatSetSizes(user->Ae, neloc, nloc, user->ne, user->n);
213:     MatSetFromOptions(user->Ae);
214:     MatSetUp(user->Ae);
215:   }

217:   MatCreate(PETSC_COMM_WORLD, &user->Ai);
218:   MatSetSizes(user->Ai, niloc, nloc, user->ni, user->n);
219:   MatSetFromOptions(user->Ai);
220:   MatSetUp(user->Ai);

222:   MatCreate(PETSC_COMM_WORLD, &user->H);
223:   MatSetSizes(user->H, nloc, nloc, user->n, user->n);
224:   MatSetFromOptions(user->H);
225:   MatSetUp(user->H);
226:   return 0;
227: }

229: PetscErrorCode DestroyProblem(AppCtx *user)
230: {
231:   if (!user->noeqflag) MatDestroy(&user->Ae);
232:   MatDestroy(&user->Ai);
233:   MatDestroy(&user->H);

235:   VecDestroy(&user->x);
236:   if (!user->noeqflag) VecDestroy(&user->ce);
237:   VecDestroy(&user->ci);
238:   VecDestroy(&user->xl);
239:   VecDestroy(&user->xu);
240:   VecDestroy(&user->Xseq);
241:   VecScatterDestroy(&user->scat);
242:   return 0;
243: }

245: /* Evaluate
246:    f(x) = (x0 - 2)^2 + (x1 - 2)^2 - 2*(x0 + x1)
247:    G = grad f = [2*(x0 - 2) - 2;
248:                  2*(x1 - 2) - 2]
249: */
250: PetscErrorCode FormFunctionGradient(Tao tao, Vec X, PetscReal *f, Vec G, void *ctx)
251: {
252:   PetscScalar        g;
253:   const PetscScalar *x;
254:   MPI_Comm           comm;
255:   PetscMPIInt        rank;
256:   PetscReal          fin;
257:   AppCtx            *user = (AppCtx *)ctx;
258:   Vec                Xseq = user->Xseq;
259:   VecScatter         scat = user->scat;

261:   PetscObjectGetComm((PetscObject)tao, &comm);
262:   MPI_Comm_rank(comm, &rank);

264:   VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD);
265:   VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD);

267:   fin = 0.0;
268:   if (rank == 0) {
269:     VecGetArrayRead(Xseq, &x);
270:     fin = (x[0] - 2.0) * (x[0] - 2.0) + (x[1] - 2.0) * (x[1] - 2.0) - 2.0 * (x[0] + x[1]);
271:     g   = 2.0 * (x[0] - 2.0) - 2.0;
272:     VecSetValue(G, 0, g, INSERT_VALUES);
273:     g = 2.0 * (x[1] - 2.0) - 2.0;
274:     VecSetValue(G, 1, g, INSERT_VALUES);
275:     VecRestoreArrayRead(Xseq, &x);
276:   }
277:   MPI_Allreduce(&fin, f, 1, MPIU_REAL, MPIU_SUM, comm);
278:   VecAssemblyBegin(G);
279:   VecAssemblyEnd(G);
280:   return 0;
281: }

283: /* Evaluate
284:    H = fxx + grad (grad g^T*DI) - grad (grad h^T*DE)]
285:      = [ 2*(1+de[0]-di[0]+di[1]), 0;
286:                    0,             2]
287: */
288: PetscErrorCode FormHessian(Tao tao, Vec x, Mat H, Mat Hpre, void *ctx)
289: {
290:   AppCtx            *user = (AppCtx *)ctx;
291:   Vec                DE, DI;
292:   const PetscScalar *de, *di;
293:   PetscInt           zero = 0, one = 1;
294:   PetscScalar        two = 2.0;
295:   PetscScalar        val = 0.0;
296:   Vec                Deseq, Diseq;
297:   VecScatter         Descat, Discat;
298:   PetscMPIInt        rank;
299:   MPI_Comm           comm;

301:   TaoGetDualVariables(tao, &DE, &DI);

303:   PetscObjectGetComm((PetscObject)tao, &comm);
304:   MPI_Comm_rank(comm, &rank);

306:   if (!user->noeqflag) {
307:     VecScatterCreateToZero(DE, &Descat, &Deseq);
308:     VecScatterBegin(Descat, DE, Deseq, INSERT_VALUES, SCATTER_FORWARD);
309:     VecScatterEnd(Descat, DE, Deseq, INSERT_VALUES, SCATTER_FORWARD);
310:   }
311:   VecScatterCreateToZero(DI, &Discat, &Diseq);
312:   VecScatterBegin(Discat, DI, Diseq, INSERT_VALUES, SCATTER_FORWARD);
313:   VecScatterEnd(Discat, DI, Diseq, INSERT_VALUES, SCATTER_FORWARD);

315:   if (rank == 0) {
316:     if (!user->noeqflag) { VecGetArrayRead(Deseq, &de); /* places equality constraint dual into array */ }
317:     VecGetArrayRead(Diseq, &di); /* places inequality constraint dual into array */

319:     if (!user->noeqflag) {
320:       val = 2.0 * (1 + de[0] - di[0] + di[1]);
321:       VecRestoreArrayRead(Deseq, &de);
322:       VecRestoreArrayRead(Diseq, &di);
323:     } else {
324:       val = 2.0 * (1 - di[0] + di[1]);
325:     }
326:     VecRestoreArrayRead(Diseq, &di);
327:     MatSetValues(H, 1, &zero, 1, &zero, &val, INSERT_VALUES);
328:     MatSetValues(H, 1, &one, 1, &one, &two, INSERT_VALUES);
329:   }
330:   MatAssemblyBegin(H, MAT_FINAL_ASSEMBLY);
331:   MatAssemblyEnd(H, MAT_FINAL_ASSEMBLY);
332:   if (!user->noeqflag) {
333:     VecScatterDestroy(&Descat);
334:     VecDestroy(&Deseq);
335:   }
336:   VecScatterDestroy(&Discat);
337:   VecDestroy(&Diseq);
338:   return 0;
339: }

341: /* Evaluate
342:    h = [ x0^2 - x1;
343:          1 -(x0^2 - x1)]
344: */
345: PetscErrorCode FormInequalityConstraints(Tao tao, Vec X, Vec CI, void *ctx)
346: {
347:   const PetscScalar *x;
348:   PetscScalar        ci;
349:   MPI_Comm           comm;
350:   PetscMPIInt        rank;
351:   AppCtx            *user = (AppCtx *)ctx;
352:   Vec                Xseq = user->Xseq;
353:   VecScatter         scat = user->scat;

355:   PetscObjectGetComm((PetscObject)tao, &comm);
356:   MPI_Comm_rank(comm, &rank);

358:   VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD);
359:   VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD);

361:   if (rank == 0) {
362:     VecGetArrayRead(Xseq, &x);
363:     ci = x[0] * x[0] - x[1];
364:     VecSetValue(CI, 0, ci, INSERT_VALUES);
365:     ci = -x[0] * x[0] + x[1] + 1.0;
366:     VecSetValue(CI, 1, ci, INSERT_VALUES);
367:     VecRestoreArrayRead(Xseq, &x);
368:   }
369:   VecAssemblyBegin(CI);
370:   VecAssemblyEnd(CI);
371:   return 0;
372: }

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;

387:   PetscObjectGetComm((PetscObject)tao, &comm);
388:   MPI_Comm_rank(comm, &rank);

390:   VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD);
391:   VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD);

393:   if (rank == 0) {
394:     VecGetArrayRead(Xseq, &x);
395:     ce = x[0] * x[0] + x[1] - 2.0;
396:     VecSetValue(CE, 0, ce, INSERT_VALUES);
397:     VecRestoreArrayRead(Xseq, &x);
398:   }
399:   VecAssemblyBegin(CE);
400:   VecAssemblyEnd(CE);
401:   return 0;
402: }

404: /*
405:   grad h = [  2*x0, -1;
406:              -2*x0,  1]
407: */
408: PetscErrorCode FormInequalityJacobian(Tao tao, Vec X, Mat JI, Mat JIpre, void *ctx)
409: {
410:   AppCtx            *user = (AppCtx *)ctx;
411:   PetscInt           zero = 0, one = 1, cols[2];
412:   PetscScalar        vals[2];
413:   const PetscScalar *x;
414:   Vec                Xseq = user->Xseq;
415:   VecScatter         scat = user->scat;
416:   MPI_Comm           comm;
417:   PetscMPIInt        rank;

419:   PetscObjectGetComm((PetscObject)tao, &comm);
420:   MPI_Comm_rank(comm, &rank);
421:   VecScatterBegin(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD);
422:   VecScatterEnd(scat, X, Xseq, INSERT_VALUES, SCATTER_FORWARD);

424:   VecGetArrayRead(Xseq, &x);
425:   if (rank == 0) {
426:     cols[0] = 0;
427:     cols[1] = 1;
428:     vals[0] = 2 * x[0];
429:     vals[1] = -1.0;
430:     MatSetValues(JI, 1, &zero, 2, cols, vals, INSERT_VALUES);
431:     vals[0] = -2 * x[0];
432:     vals[1] = 1.0;
433:     MatSetValues(JI, 1, &one, 2, cols, vals, INSERT_VALUES);
434:   }
435:   VecRestoreArrayRead(Xseq, &x);
436:   MatAssemblyBegin(JI, MAT_FINAL_ASSEMBLY);
437:   MatAssemblyEnd(JI, MAT_FINAL_ASSEMBLY);
438:   return 0;
439: }

441: /*
442:   grad g = [2*x0
443:              1.0 ]
444: */
445: PetscErrorCode FormEqualityJacobian(Tao tao, Vec X, Mat JE, Mat JEpre, void *ctx)
446: {
447:   PetscInt           zero = 0, cols[2];
448:   PetscScalar        vals[2];
449:   const PetscScalar *x;
450:   PetscMPIInt        rank;
451:   MPI_Comm           comm;

453:   PetscObjectGetComm((PetscObject)tao, &comm);
454:   MPI_Comm_rank(comm, &rank);

456:   if (rank == 0) {
457:     VecGetArrayRead(X, &x);
458:     cols[0] = 0;
459:     cols[1] = 1;
460:     vals[0] = 2 * x[0];
461:     vals[1] = 1.0;
462:     MatSetValues(JE, 1, &zero, 2, cols, vals, INSERT_VALUES);
463:     VecRestoreArrayRead(X, &x);
464:   }
465:   MatAssemblyBegin(JE, MAT_FINAL_ASSEMBLY);
466:   MatAssemblyEnd(JE, MAT_FINAL_ASSEMBLY);
467:   return 0;
468: }

470: /*TEST

472:    build:
473:       requires: !complex !defined(PETSC_USE_CXX)

475:    test:
476:       args: -tao_converged_reason -tao_gatol 1.e-6 -tao_type pdipm -tao_pdipm_kkt_shift_pd
477:       requires: mumps
478:       filter: sed  -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"

480:    test:
481:       suffix: 2
482:       args: -tao_converged_reason
483:       filter: sed  -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"

485:    test:
486:       suffix: 3
487:       args: -tao_converged_reason -no_eq
488:       filter: sed  -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"

490:    test:
491:       suffix: 4
492:       args: -tao_converged_reason -tao_almm_type classic
493:       requires: !single
494:       filter: sed  -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"

496:    test:
497:       suffix: 5
498:       args: -tao_converged_reason -tao_almm_type classic -no_eq
499:       requires: !single
500:       filter: sed  -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"

502:    test:
503:       suffix: 6
504:       args: -tao_converged_reason -tao_almm_subsolver_tao_type bqnktr
505:       filter: sed  -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"

507:    test:
508:       suffix: 7
509:       args: -tao_converged_reason -tao_almm_subsolver_tao_type bncg
510:       filter: sed  -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"

512:    test:
513:       suffix: 8
514:       nsize: 2
515:       args: -tao_converged_reason
516:       filter: sed  -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"

518:    test:
519:       suffix: 9
520:       nsize: 2
521:       args: -tao_converged_reason -vec_type cuda -mat_type aijcusparse
522:       requires: cuda
523:       filter: sed  -e "s/CONVERGED_GATOL iterations *[0-9]\{1,\}/CONVERGED_GATOL/g"

525: TEST*/