Actual source code: ex34.c

slepc-3.16.1 2021-11-17
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2021, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    This is a nonlinear eigenvalue problem. When p=2, it is reduced to a linear Laplace eigenvalue
 12:    problem.

 14:    -\nabla\cdot(|\nabla u|^{p-2} \nabla u) = k |u|^{p-2} u in (0,1)x(0,1),

 16:    u = 0 on the entire boundary.

 18:    The code is implemented based on DMPlex using Q1 FEM on a quadrilateral mesh. In this code, we consider p=3.

 20:    Contributed  by Fande Kong fdkong.jd@gmail.com
 21: */

 23: static char help[] = "Nonlinear inverse iteration for A(x)*x=lambda*B(x)*x.\n\n";

 25: #include <slepceps.h>
 26: #include <petscdmplex.h>
 27: #include <petscds.h>

 29: PetscErrorCode CreateSquareMesh(MPI_Comm,DM*);
 30: PetscErrorCode SetupDiscretization(DM);
 31: PetscErrorCode FormJacobianA(SNES,Vec,Mat,Mat,void*);
 32: PetscErrorCode FormFunctionA(SNES,Vec,Vec,void*);
 33: PetscErrorCode MatMult_A(Mat A,Vec x,Vec y);
 34: PetscErrorCode FormJacobianB(SNES,Vec,Mat,Mat,void*);
 35: PetscErrorCode FormFunctionB(SNES,Vec,Vec,void*);
 36: PetscErrorCode MatMult_B(Mat A,Vec x,Vec y);
 37: PetscErrorCode FormFunctionAB(SNES,Vec,Vec,Vec,void*);
 38: PetscErrorCode BoundaryGlobalIndex(DM,const char*,IS*);

 40: typedef struct {
 41:   IS    bdis; /* global indices for boundary DoFs */
 42:   SNES  snes;
 43: } AppCtx;

 45: int main(int argc,char **argv)
 46: {
 47:   DM             dm;
 48:   MPI_Comm       comm;
 49:   AppCtx         user;
 50:   EPS            eps;  /* eigenproblem solver context */
 51:   ST             st;
 52:   EPSType        type;
 53:   Mat            A,B,P;
 54:   Vec            v0;
 55:   PetscContainer container;
 56:   PetscInt       nev,nconv,m,n,M,N;
 57:   PetscBool      nonlin,flg=PETSC_FALSE,update;
 58:   SNES           snes;
 59:   PetscReal      tol,relerr;
 60:   PetscBool      use_shell_matrix=PETSC_FALSE,test_init_sol=PETSC_FALSE;

 63:   SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
 64:   comm = PETSC_COMM_WORLD;
 65:   /* Create a quadrilateral mesh on domain (0,1)x(0,1) */
 66:   CreateSquareMesh(comm,&dm);
 67:   /* Setup basis function */
 68:   SetupDiscretization(dm);
 69:   BoundaryGlobalIndex(dm,"marker",&user.bdis);
 70:   /* Check if we are going to use shell matrices */
 71:   PetscOptionsGetBool(NULL,NULL,"-use_shell_matrix",&use_shell_matrix,NULL);
 72:   if (use_shell_matrix) {
 73:     DMCreateMatrix(dm,&P);
 74:     MatGetLocalSize(P,&m,&n);
 75:     MatGetSize(P,&M,&N);
 76:     MatCreateShell(comm,m,n,M,N,&user,&A);
 77:     MatShellSetOperation(A,MATOP_MULT,(void(*)(void))MatMult_A);
 78:     MatCreateShell(comm,m,n,M,N,&user,&B);
 79:     MatShellSetOperation(B,MATOP_MULT,(void(*)(void))MatMult_B);
 80:   } else {
 81:     DMCreateMatrix(dm,&A);
 82:     MatDuplicate(A,MAT_COPY_VALUES,&B);
 83:   }

 85:   /*
 86:      Compose callback functions and context that will be needed by the solver
 87:   */
 88:   PetscObjectComposeFunction((PetscObject)A,"formFunction",FormFunctionA);
 89:   PetscOptionsGetBool(NULL,NULL,"-form_function_ab",&flg,NULL);
 90:   if (flg) {
 91:     PetscObjectComposeFunction((PetscObject)A,"formFunctionAB",FormFunctionAB);
 92:   }
 93:   PetscObjectComposeFunction((PetscObject)A,"formJacobian",FormJacobianA);
 94:   PetscObjectComposeFunction((PetscObject)B,"formFunction",FormFunctionB);
 95:   PetscContainerCreate(comm,&container);
 96:   PetscContainerSetPointer(container,&user);
 97:   PetscObjectCompose((PetscObject)A,"formFunctionCtx",(PetscObject)container);
 98:   PetscObjectCompose((PetscObject)A,"formJacobianCtx",(PetscObject)container);
 99:   PetscObjectCompose((PetscObject)B,"formFunctionCtx",(PetscObject)container);
100:   PetscContainerDestroy(&container);

102:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
103:                 Create the eigensolver and set various options
104:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

106:   EPSCreate(comm,&eps);
107:   EPSSetOperators(eps,A,B);
108:   EPSSetProblemType(eps,EPS_GNHEP);
109:   /*
110:      Use nonlinear inverse iteration
111:   */
112:   EPSSetType(eps,EPSPOWER);
113:   EPSPowerSetNonlinear(eps,PETSC_TRUE);
114:   /*
115:     Attach DM to SNES
116:   */
117:   EPSPowerGetSNES(eps,&snes);
118:   user.snes = snes;
119:   SNESSetDM(snes,dm);
120:   EPSSetFromOptions(eps);

122:   /* Set a preconditioning matrix to ST */
123:   if (use_shell_matrix) {
124:     EPSGetST(eps,&st);
125:     STSetPreconditionerMat(st,P);
126:   }

128:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
129:                       Solve the eigensystem
130:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

132:   EPSSolve(eps);

134:   EPSGetConverged(eps,&nconv);
135:   PetscOptionsGetBool(NULL,NULL,"-test_init_sol",&test_init_sol,NULL);
136:   if (nconv && test_init_sol) {
137:     PetscScalar   k;
138:     PetscReal     norm0;
139:     PetscInt      nits;

141:     MatCreateVecs(A,&v0,NULL);
142:     EPSGetEigenpair(eps,0,&k,NULL,v0,NULL);
143:     EPSSetInitialSpace(eps,1,&v0);
144:     VecDestroy(&v0);
145:     /* Norm of the previous residual */
146:     SNESGetFunctionNorm(snes,&norm0);
147:     /* Make the tolerance smaller than the last residual
148:        SNES will converge right away if the initial is setup correctly */
149:     SNESSetTolerances(snes,norm0*1.2,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
150:     EPSSolve(eps);
151:     /* Number of Newton iterations supposes to be zero */
152:     SNESGetIterationNumber(snes,&nits);
153:     if (nits) {
154:       PetscPrintf(comm," Number of Newton iterations %D should be zero \n",nits);
155:     }
156:   }

158:   /*
159:      Optional: Get some information from the solver and display it
160:   */
161:   EPSGetType(eps,&type);
162:   EPSGetTolerances(eps,&tol,NULL);
163:   EPSPowerGetNonlinear(eps,&nonlin);
164:   EPSPowerGetUpdate(eps,&update);
165:   PetscPrintf(comm," Solution method: %s%s\n\n",type,nonlin?(update?" (nonlinear with monolithic update)":" (nonlinear)"):"");
166:   EPSGetDimensions(eps,&nev,NULL,NULL);
167:   PetscPrintf(comm," Number of requested eigenvalues: %D\n",nev);

169:   /* print eigenvalue and error */
170:   EPSGetConverged(eps,&nconv);
171:   if (nconv>0) {
172:     PetscScalar   k;
173:     PetscReal     na,nb;
174:     Vec           a,b,eigen;
175:     DMCreateGlobalVector(dm,&a);
176:     VecDuplicate(a,&b);
177:     VecDuplicate(a,&eigen);
178:     EPSGetEigenpair(eps,0,&k,NULL,eigen,NULL);
179:     FormFunctionA(snes,eigen,a,&user);
180:     FormFunctionB(snes,eigen,b,&user);
181:     VecAXPY(a,-k,b);
182:     VecNorm(a,NORM_2,&na);
183:     VecNorm(b,NORM_2,&nb);
184:     relerr = na/(nb*PetscAbsScalar(k));
185:     if (relerr<10*tol) {
186:       PetscPrintf(comm,"k: %g, relative error below tol\n",(double)PetscRealPart(k));
187:     } else {
188:       PetscPrintf(comm,"k: %g, relative error: %g\n",(double)PetscRealPart(k),(double)relerr);
189:     }
190:     VecDestroy(&a);
191:     VecDestroy(&b);
192:     VecDestroy(&eigen);
193:   } else {
194:     PetscPrintf(comm,"Solver did not converge\n");
195:   }

197:   MatDestroy(&A);
198:   MatDestroy(&B);
199:   if (use_shell_matrix) {
200:     MatDestroy(&P);
201:   }
202:   DMDestroy(&dm);
203:   EPSDestroy(&eps);
204:   ISDestroy(&user.bdis);
205:   SlepcFinalize();
206:   return ierr;
207: }

209: /* <|u|u, v> */
210: static void f0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
211:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
212:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
213:                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
214: {
215:   PetscScalar cof = PetscAbsScalar(u[0]);

217:   f0[0] = cof*u[0];
218: }

220: /* <|\nabla u| \nabla u, \nabla v> */
221: static void f1_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
222:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
223:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
224:                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
225: {
226:   PetscInt    d;
227:   PetscScalar cof = 0;
228:   for (d = 0; d < dim; ++d)  cof += u_x[d]*u_x[d];

230:   cof = PetscSqrtScalar(cof);

232:   for (d = 0; d < dim; ++d) f1[d] = u_x[d]*cof;
233: }

235: /* approximate  Jacobian for   <|u|u, v> */
236: static void g0_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
237:                   const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
238:                   const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
239:                   PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[])
240: {
241:   g0[0] = 1.0*PetscAbsScalar(u[0]);
242: }

244: /* approximate  Jacobian for   <|\nabla u| \nabla u, \nabla v> */
245: static void g3_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
246:                   const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
247:                   const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
248:                   PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
249: {
250:   PetscInt d;

252:   for (d = 0; d < dim; ++d) g3[d*dim+d] = 1.0;
253: }

255: PetscErrorCode SetupDiscretization(DM dm)
256: {
257:   PetscFE        fe;
258:   MPI_Comm       comm;

262:   /* Create finite element */
263:   PetscObjectGetComm((PetscObject)dm,&comm);
264:   PetscFECreateDefault(comm,2,1,PETSC_FALSE,NULL,-1,&fe);
265:   PetscObjectSetName((PetscObject)fe,"u");
266:   DMSetField(dm,0,NULL,(PetscObject)fe);
267:   DMCreateDS(dm);
268:   PetscFEDestroy(&fe);
269:   return(0);
270: }

272: PetscErrorCode CreateSquareMesh(MPI_Comm comm,DM *dm)
273: {
274:   PetscInt       cells[] = {8,8};
275:   PetscInt       dim = 2;
276:   DM             pdm;
277:   PetscMPIInt    size;

281:   DMPlexCreateBoxMesh(comm,dim,PETSC_FALSE,cells,NULL,NULL,NULL,PETSC_TRUE,dm);
282:   DMSetFromOptions(*dm);
283:   DMSetUp(*dm);
284:   MPI_Comm_size(comm,&size);
285:   if (size > 1) {
286:     DMPlexDistribute(*dm,0,NULL,&pdm);
287:     DMDestroy(dm);
288:     *dm = pdm;
289:   }
290:   return(0);
291: }

293: PetscErrorCode BoundaryGlobalIndex(DM dm,const char labelname[],IS *bdis)
294: {
295:   IS             bdpoints;
296:   PetscInt       nindices,*indices,numDof,offset,npoints,i,j;
297:   const PetscInt *bdpoints_indices;
298:   DMLabel        bdmarker;
299:   PetscSection   gsection;

303:   DMGetGlobalSection(dm,&gsection);
304:   DMGetLabel(dm,labelname,&bdmarker);
305:   DMLabelGetStratumIS(bdmarker,1,&bdpoints);
306:   ISGetLocalSize(bdpoints,&npoints);
307:   ISGetIndices(bdpoints,&bdpoints_indices);
308:   nindices = 0;
309:   for (i=0;i<npoints;i++) {
310:     PetscSectionGetDof(gsection,bdpoints_indices[i],&numDof);
311:     if (numDof<=0) continue;
312:     nindices += numDof;
313:   }
314:   PetscCalloc1(nindices,&indices);
315:   nindices = 0;
316:   for (i=0;i<npoints;i++) {
317:     PetscSectionGetDof(gsection,bdpoints_indices[i],&numDof);
318:     if (numDof<=0) continue;
319:     PetscSectionGetOffset(gsection,bdpoints_indices[i],&offset);
320:     for (j=0;j<numDof;j++) indices[nindices++] = offset+j;
321:   }
322:   ISRestoreIndices(bdpoints,&bdpoints_indices);
323:   ISDestroy(&bdpoints);
324:   ISCreateGeneral(PetscObjectComm((PetscObject)dm),nindices,indices,PETSC_OWN_POINTER,bdis);
325:   return(0);
326: }

328: static PetscErrorCode FormJacobian(SNES snes,Vec X,Mat A,Mat B,void *ctx)
329: {
330:   DM             dm;
331:   Vec            Xloc;

335:   SNESGetDM(snes,&dm);
336:   DMGetLocalVector(dm,&Xloc);
337:   VecZeroEntries(Xloc);
338:   DMGlobalToLocalBegin(dm,X,INSERT_VALUES,Xloc);
339:   DMGlobalToLocalEnd(dm,X,INSERT_VALUES,Xloc);
340:   CHKMEMQ;
341:   DMPlexSNESComputeJacobianFEM(dm,Xloc,A,B,ctx);
342:   if (A!=B) {
343:     MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
344:     MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
345:   }
346:   CHKMEMQ;
347:   DMRestoreLocalVector(dm,&Xloc);
348:   return(0);
349: }

351: PetscErrorCode FormJacobianA(SNES snes,Vec X,Mat A,Mat B,void *ctx)
352: {
354:   DM             dm;
355:   PetscDS        prob;
356:   PetscWeakForm  wf;
357:   AppCtx         *userctx = (AppCtx *)ctx;

360:   MatSetOption(B,MAT_KEEP_NONZERO_PATTERN,PETSC_TRUE);
361:   SNESGetDM(snes,&dm);
362:   DMGetDS(dm,&prob);
363:   PetscDSGetWeakForm(prob, &wf);
364:   PetscWeakFormClearIndex(wf, NULL, 0, 0, 0, PETSC_WF_G3, 0);
365:   PetscWeakFormSetIndexJacobian(wf, NULL, 0, 0, 0, 0, 0, NULL, 0, NULL, 0, NULL, 0, g3_uu);
366:   FormJacobian(snes,X,A,B,ctx);
367:   MatZeroRowsIS(B,userctx->bdis,1.0,NULL,NULL);
368:   return(0);
369: }

371: PetscErrorCode FormJacobianB(SNES snes,Vec X,Mat A,Mat B,void *ctx)
372: {
374:   DM             dm;
375:   PetscDS        prob;
376:   PetscWeakForm  wf;
377:   AppCtx         *userctx = (AppCtx *)ctx;

380:   MatSetOption(B,MAT_KEEP_NONZERO_PATTERN,PETSC_TRUE);
381:   SNESGetDM(snes,&dm);
382:   DMGetDS(dm,&prob);
383:   PetscDSGetWeakForm(prob, &wf);
384:   PetscWeakFormClearIndex(wf, NULL, 0, 0, 0, PETSC_WF_G3, 0);
385:   PetscWeakFormSetIndexJacobian(wf, NULL, 0, 0, 0, 0, 0, g0_uu, 0, NULL, 0, NULL, 0, NULL);
386:   FormJacobian(snes,X,A,B,ctx);
387:   MatZeroRowsIS(B,userctx->bdis,0.0,NULL,NULL);
388:   return(0);
389: }

391: PetscErrorCode FormFunctionAB(SNES snes,Vec x,Vec Ax,Vec Bx,void *ctx)
392: {

396:   /*
397:    * In real applications, users should have a generic formFunctionAB which
398:    * forms Ax and Bx simultaneously for an more efficient calculation.
399:    * In this example, we just call FormFunctionA+FormFunctionB to mimic how
400:    * to use FormFunctionAB
401:    */
402:   FormFunctionA(snes,x,Ax,ctx);
403:   FormFunctionB(snes,x,Bx,ctx);
404:   return(0);
405: }

407: static PetscErrorCode FormFunction(SNES snes,Vec X,Vec F,void *ctx)
408: {
409:   DM             dm;
410:   Vec            Xloc,Floc;

414:   SNESGetDM(snes,&dm);
415:   DMGetLocalVector(dm,&Xloc);
416:   DMGetLocalVector(dm,&Floc);
417:   VecZeroEntries(Xloc);
418:   VecZeroEntries(Floc);
419:   DMGlobalToLocalBegin(dm,X,INSERT_VALUES,Xloc);
420:   DMGlobalToLocalEnd(dm,X,INSERT_VALUES,Xloc);
421:   CHKMEMQ;
422:   DMPlexSNESComputeResidualFEM(dm,Xloc,Floc,ctx);
423:   CHKMEMQ;
424:   VecZeroEntries(F);
425:   DMLocalToGlobalBegin(dm,Floc,ADD_VALUES,F);
426:   DMLocalToGlobalEnd(dm,Floc,ADD_VALUES,F);
427:   DMRestoreLocalVector(dm,&Xloc);
428:   DMRestoreLocalVector(dm,&Floc);
429:   return(0);
430: }

432: PetscErrorCode FormFunctionA(SNES snes,Vec X,Vec F,void *ctx)
433: {
435:   DM             dm;
436:   PetscDS        prob;
437:   PetscWeakForm  wf;
438:   PetscInt       nindices,iStart,iEnd,i;
439:   AppCtx         *userctx = (AppCtx *)ctx;
440:   PetscScalar    *array,value;
441:   const PetscInt *indices;
442:   PetscInt       vecstate;

445:   SNESGetDM(snes,&dm);
446:   DMGetDS(dm,&prob);
447:   /* hook functions */
448:   PetscDSGetWeakForm(prob, &wf);
449:   PetscWeakFormClearIndex(wf, NULL, 0, 0, 0, PETSC_WF_F0, 0);
450:   PetscWeakFormSetIndexResidual(wf, NULL, 0, 0, 0, 0, NULL, 0, f1_u);
451:   FormFunction(snes,X,F,ctx);
452:   /* Boundary condition */
453:   VecLockGet(X,&vecstate);
454:   if (vecstate>0) {
455:     VecLockReadPop(X);
456:   }
457:   VecGetOwnershipRange(X,&iStart,&iEnd);
458:   VecGetArray(X,&array);
459:   ISGetLocalSize(userctx->bdis,&nindices);
460:   ISGetIndices(userctx->bdis,&indices);
461:   for (i=0;i<nindices;i++) {
462:     value = array[indices[i]-iStart] - 0.0;
463:     VecSetValue(F,indices[i],value,INSERT_VALUES);
464:   }
465:   ISRestoreIndices(userctx->bdis,&indices);
466:   VecRestoreArray(X,&array);
467:   if (vecstate>0) {
468:     VecLockReadPush(X);
469:   }
470:   VecAssemblyBegin(F);
471:   VecAssemblyEnd(F);
472:   return(0);
473: }

475: PetscErrorCode MatMult_A(Mat A,Vec x,Vec y)
476: {
478:   AppCtx         *userctx;

481:   MatShellGetContext(A,&userctx);
482:   FormFunctionA(userctx->snes,x,y,userctx);
483:   return(0);
484: }

486: PetscErrorCode FormFunctionB(SNES snes,Vec X,Vec F,void *ctx)
487: {
489:   DM             dm;
490:   PetscDS        prob;
491:   PetscWeakForm  wf;
492:   PetscInt       nindices,iStart,iEnd,i;
493:   AppCtx         *userctx = (AppCtx *)ctx;
494:   PetscScalar    value;
495:   const PetscInt *indices;

498:   SNESGetDM(snes,&dm);
499:   DMGetDS(dm,&prob);
500:   /* hook functions */
501:   PetscDSGetWeakForm(prob, &wf);
502:   PetscWeakFormClearIndex(wf, NULL, 0, 0, 0, PETSC_WF_F1, 0);
503:   PetscWeakFormSetIndexResidual(wf, NULL, 0, 0, 0, 0, f0_u, 0, NULL);
504:   FormFunction(snes,X,F,ctx);
505:   /* Boundary condition */
506:   VecGetOwnershipRange(F,&iStart,&iEnd);
507:   ISGetLocalSize(userctx->bdis,&nindices);
508:   ISGetIndices(userctx->bdis,&indices);
509:   for (i=0;i<nindices;i++) {
510:     value = 0.0;
511:     VecSetValue(F,indices[i],value,INSERT_VALUES);
512:   }
513:   ISRestoreIndices(userctx->bdis,&indices);
514:   VecAssemblyBegin(F);
515:   VecAssemblyEnd(F);
516:   return(0);
517: }

519: PetscErrorCode MatMult_B(Mat B,Vec x,Vec y)
520: {
522:   AppCtx         *userctx;

525:   MatShellGetContext(B,&userctx);
526:   FormFunctionB(userctx->snes,x,y,userctx);
527:   return(0);
528: }

530: /*TEST

532:    testset:
533:       requires: double
534:       args: -petscspace_degree 1 -petscspace_poly_tensor
535:       output_file: output/ex34_1.out
536:       test:
537:          suffix: 1
538:       test:
539:          suffix: 2
540:          args: -eps_power_update -form_function_ab {{0 1}}
541:          filter: sed -e "s/ with monolithic update//"
542:       test:
543:          suffix: 3
544:          args: -use_shell_matrix -eps_power_snes_mf_operator 1
545:       test:
546:          suffix: 4
547:          args: -use_shell_matrix -eps_power_update -init_eps_power_snes_mf_operator 1 -eps_power_snes_mf_operator 1 -form_function_ab {{0 1}}
548:          filter: sed -e "s/ with monolithic update//"
549:       test:
550:          suffix: 5
551:          args: -use_shell_matrix -eps_power_update -init_eps_power_snes_mf_operator 1 -eps_power_snes_mf_operator 1 -form_function_ab {{0 1}} -test_init_sol 1
552:          filter: sed -e "s/ with monolithic update//"

554:       test:
555:          suffix: 6
556:          args: -use_shell_matrix -eps_power_update -init_eps_power_snes_mf_operator 1 -eps_power_snes_mf_operator 1 -form_function_ab {{0 1}} -eps_monitor_all
557:          output_file: output/ex34_6.out
558:          filter: sed -e "s/\([+-].*i\)//g" -e "1,3s/[0-9]//g" -e "/[45] EPS/d"
559: TEST*/