Actual source code: ex5.c

  1: static char help[] = "Bratu nonlinear PDE in 2d.\n\
  2: We solve the  Bratu (SFI - solid fuel ignition) problem in a 2D rectangular\n\
  3: domain, using distributed arrays (DMDAs) to partition the parallel grid.\n\
  4: The command line options include:\n\
  5:   -par <parameter>, where <parameter> indicates the problem's nonlinearity\n\
  6:      problem SFI:  <parameter> = Bratu parameter (0 <= par <= 6.81)\n\n\
  7:   -m_par/n_par <parameter>, where <parameter> indicates an integer\n \
  8:       that MMS3 will be evaluated with 2^m_par, 2^n_par";

 10: /*T
 11:    Concepts: SNES^parallel Bratu example
 12:    Concepts: DMDA^using distributed arrays;
 13:    Concepts: IS coloirng types;
 14:    Processors: n
 15: T*/

 17: /* ------------------------------------------------------------------------

 19:     Solid Fuel Ignition (SFI) problem.  This problem is modeled by
 20:     the partial differential equation

 22:             -Laplacian u - lambda*exp(u) = 0,  0 < x,y < 1,

 24:     with boundary conditions

 26:              u = 0  for  x = 0, x = 1, y = 0, y = 1.

 28:     A finite difference approximation with the usual 5-point stencil
 29:     is used to discretize the boundary value problem to obtain a nonlinear
 30:     system of equations.

 32:       This example shows how geometric multigrid can be run transparently with a nonlinear solver so long
 33:       as SNESSetDM() is provided. Example usage

 35:       ./ex5 -pc_type mg -ksp_monitor  -snes_view -pc_mg_levels 3 -pc_mg_galerkin pmat -da_grid_x 17 -da_grid_y 17
 36:              -mg_levels_ksp_monitor -snes_monitor -mg_levels_pc_type sor -pc_mg_type full

 38:       or to run with grid sequencing on the nonlinear problem (note that you do not need to provide the number of
 39:          multigrid levels, it will be determined automatically based on the number of refinements done)

 41:       ./ex5 -pc_type mg -ksp_monitor  -snes_view -pc_mg_galerkin pmat -snes_grid_sequence 3
 42:              -mg_levels_ksp_monitor -snes_monitor -mg_levels_pc_type sor -pc_mg_type full

 44:   ------------------------------------------------------------------------- */

 46: /*
 47:    Include "petscdmda.h" so that we can use distributed arrays (DMDAs).
 48:    Include "petscsnes.h" so that we can use SNES solvers.  Note that this
 49: */
 50: #include <petscdm.h>
 51: #include <petscdmda.h>
 52: #include <petscsnes.h>
 53: #include <petscmatlab.h>
 54: #include <petsc/private/snesimpl.h>

 56: /*
 57:    User-defined application context - contains data needed by the
 58:    application-provided call-back routines, FormJacobianLocal() and
 59:    FormFunctionLocal().
 60: */
 61: typedef struct AppCtx AppCtx;
 62: struct AppCtx {
 63:   PetscReal param;          /* test problem parameter */
 64:   PetscInt  m,n;            /* MMS3 parameters */
 65:   PetscErrorCode (*mms_solution)(AppCtx*,const DMDACoor2d*,PetscScalar*);
 66:   PetscErrorCode (*mms_forcing)(AppCtx*,const DMDACoor2d*,PetscScalar*);
 67: };

 69: /*
 70:    User-defined routines
 71: */
 72: extern PetscErrorCode FormInitialGuess(DM,AppCtx*,Vec);
 73: extern PetscErrorCode FormFunctionLocal(DMDALocalInfo*,PetscScalar**,PetscScalar**,AppCtx*);
 74: extern PetscErrorCode FormExactSolution(DM,AppCtx*,Vec);
 75: extern PetscErrorCode ZeroBCSolution(AppCtx*,const DMDACoor2d*,PetscScalar*);
 76: extern PetscErrorCode MMSSolution1(AppCtx*,const DMDACoor2d*,PetscScalar*);
 77: extern PetscErrorCode MMSForcing1(AppCtx*,const DMDACoor2d*,PetscScalar*);
 78: extern PetscErrorCode MMSSolution2(AppCtx*,const DMDACoor2d*,PetscScalar*);
 79: extern PetscErrorCode MMSForcing2(AppCtx*,const DMDACoor2d*,PetscScalar*);
 80: extern PetscErrorCode MMSSolution3(AppCtx*,const DMDACoor2d*,PetscScalar*);
 81: extern PetscErrorCode MMSForcing3(AppCtx*,const DMDACoor2d*,PetscScalar*);
 82: extern PetscErrorCode MMSSolution4(AppCtx*,const DMDACoor2d*,PetscScalar*);
 83: extern PetscErrorCode MMSForcing4(AppCtx*,const DMDACoor2d*,PetscScalar*);
 84: extern PetscErrorCode FormJacobianLocal(DMDALocalInfo*,PetscScalar**,Mat,Mat,AppCtx*);
 85: extern PetscErrorCode FormObjectiveLocal(DMDALocalInfo*,PetscScalar**,PetscReal*,AppCtx*);
 86: extern PetscErrorCode FormFunctionMatlab(SNES,Vec,Vec,void*);
 87: extern PetscErrorCode NonlinearGS(SNES,Vec,Vec,void*);

 89: int main(int argc,char **argv)
 90: {
 91:   SNES           snes;                         /* nonlinear solver */
 92:   Vec            x;                            /* solution vector */
 93:   AppCtx         user;                         /* user-defined work context */
 94:   PetscInt       its;                          /* iterations for convergence */
 96:   PetscReal      bratu_lambda_max = 6.81;
 97:   PetscReal      bratu_lambda_min = 0.;
 98:   PetscInt       MMS              = 0;
 99:   PetscBool      flg              = PETSC_FALSE;
100:   DM             da;
101:   Vec            r               = NULL;
102:   KSP            ksp;
103:   PetscInt       lits,slits;

105:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
106:      Initialize program
107:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

109:   PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;

111:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
112:      Initialize problem parameters
113:   - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
114:   user.param = 6.0;
115:   PetscOptionsGetReal(NULL,NULL,"-par",&user.param,NULL);
116:   if (user.param > bratu_lambda_max || user.param < bratu_lambda_min) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lambda, %g, is out of range, [%g, %g]", user.param, bratu_lambda_min, bratu_lambda_max);
117:   PetscOptionsGetInt(NULL,NULL,"-mms",&MMS,NULL);
118:   if (MMS == 3) {
119:     PetscInt mPar = 2, nPar = 1;
120:     PetscOptionsGetInt(NULL,NULL,"-m_par",&mPar,NULL);
121:     PetscOptionsGetInt(NULL,NULL,"-n_par",&nPar,NULL);
122:     user.m = PetscPowInt(2,mPar);
123:     user.n = PetscPowInt(2,nPar);
124:   }

126:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
127:      Create nonlinear solver context
128:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
129:   SNESCreate(PETSC_COMM_WORLD,&snes);
130:   SNESSetCountersReset(snes,PETSC_FALSE);
131:   SNESSetNGS(snes, NonlinearGS, NULL);

133:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
134:      Create distributed array (DMDA) to manage parallel grid and vectors
135:   - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
136:   DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE,DMDA_STENCIL_STAR,4,4,PETSC_DECIDE,PETSC_DECIDE,1,1,NULL,NULL,&da);
137:   DMSetFromOptions(da);
138:   DMSetUp(da);
139:   DMDASetUniformCoordinates(da, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0);
140:   DMSetApplicationContext(da,&user);
141:   SNESSetDM(snes,da);
142:   /*  - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
143:      Extract global vectors from DMDA; then duplicate for remaining
144:      vectors that are the same types
145:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
146:   DMCreateGlobalVector(da,&x);

148:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
149:      Set local function evaluation routine
150:   - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
151:   user.mms_solution = ZeroBCSolution;
152:   switch (MMS) {
153:   case 0: user.mms_solution = NULL; user.mms_forcing = NULL;
154:   case 1: user.mms_solution = MMSSolution1; user.mms_forcing = MMSForcing1; break;
155:   case 2: user.mms_solution = MMSSolution2; user.mms_forcing = MMSForcing2; break;
156:   case 3: user.mms_solution = MMSSolution3; user.mms_forcing = MMSForcing3; break;
157:   case 4: user.mms_solution = MMSSolution4; user.mms_forcing = MMSForcing4; break;
158:   default: SETERRQ1(PETSC_COMM_WORLD,PETSC_ERR_USER,"Unknown MMS type %d",MMS);
159:   }
160:   DMDASNESSetFunctionLocal(da,INSERT_VALUES,(DMDASNESFunction)FormFunctionLocal,&user);
161:   PetscOptionsGetBool(NULL,NULL,"-fd",&flg,NULL);
162:   if (!flg) {
163:     DMDASNESSetJacobianLocal(da,(DMDASNESJacobian)FormJacobianLocal,&user);
164:   }

166:   PetscOptionsGetBool(NULL,NULL,"-obj",&flg,NULL);
167:   if (flg) {
168:     DMDASNESSetObjectiveLocal(da,(DMDASNESObjective)FormObjectiveLocal,&user);
169:   }

171:   if (PetscDefined(HAVE_MATLAB_ENGINE)) {
172:     PetscBool matlab_function = PETSC_FALSE;
173:     PetscOptionsGetBool(NULL,NULL,"-matlab_function",&matlab_function,0);
174:     if (matlab_function) {
175:       VecDuplicate(x,&r);
176:       SNESSetFunction(snes,r,FormFunctionMatlab,&user);
177:     }
178:   }

180:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
181:      Customize nonlinear solver; set runtime options
182:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
183:   SNESSetFromOptions(snes);

185:   FormInitialGuess(da,&user,x);

187:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
188:      Solve nonlinear system
189:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
190:   SNESSolve(snes,NULL,x);
191:   SNESGetIterationNumber(snes,&its);

193:   SNESGetLinearSolveIterations(snes,&slits);
194:   SNESGetKSP(snes,&ksp);
195:   KSPGetTotalIterations(ksp,&lits);
196:   if (lits != slits) SETERRQ2(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Number of total linear iterations reported by SNES %D does not match reported by KSP %D",slits,lits);
197:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
198:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

200:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
201:      If using MMS, check the l_2 error
202:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
203:   if (MMS) {
204:     Vec       e;
205:     PetscReal errorl2, errorinf;
206:     PetscInt  N;

208:     VecDuplicate(x, &e);
209:     PetscObjectViewFromOptions((PetscObject) x, NULL, "-sol_view");
210:     FormExactSolution(da, &user, e);
211:     PetscObjectViewFromOptions((PetscObject) e, NULL, "-exact_view");
212:     VecAXPY(e, -1.0, x);
213:     PetscObjectViewFromOptions((PetscObject) e, NULL, "-error_view");
214:     VecNorm(e, NORM_2, &errorl2);
215:     VecNorm(e, NORM_INFINITY, &errorinf);
216:     VecGetSize(e, &N);
217:     PetscPrintf(PETSC_COMM_WORLD, "N: %D error L2 %g inf %g\n", N, (double) errorl2/PetscSqrtReal(N), (double) errorinf);
218:     VecDestroy(&e);
219:     PetscLogEventSetDof(SNES_Solve, 0, N);
220:     PetscLogEventSetError(SNES_Solve, 0, errorl2/PetscSqrtReal(N));
221:   }

223:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
224:      Free work space.  All PETSc objects should be destroyed when they
225:      are no longer needed.
226:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
227:   VecDestroy(&r);
228:   VecDestroy(&x);
229:   SNESDestroy(&snes);
230:   DMDestroy(&da);
231:   PetscFinalize();
232:   return ierr;
233: }
234: /* ------------------------------------------------------------------- */
235: /*
236:    FormInitialGuess - Forms initial approximation.

238:    Input Parameters:
239:    da - The DM
240:    user - user-defined application context

242:    Output Parameter:
243:    X - vector
244:  */
245: PetscErrorCode FormInitialGuess(DM da,AppCtx *user,Vec X)
246: {
247:   PetscInt       i,j,Mx,My,xs,ys,xm,ym;
249:   PetscReal      lambda,temp1,temp,hx,hy;
250:   PetscScalar    **x;

253:   DMDAGetInfo(da,PETSC_IGNORE,&Mx,&My,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE);

255:   lambda = user->param;
256:   hx     = 1.0/(PetscReal)(Mx-1);
257:   hy     = 1.0/(PetscReal)(My-1);
258:   temp1  = lambda/(lambda + 1.0);

260:   /*
261:      Get a pointer to vector data.
262:        - For default PETSc vectors, VecGetArray() returns a pointer to
263:          the data array.  Otherwise, the routine is implementation dependent.
264:        - You MUST call VecRestoreArray() when you no longer need access to
265:          the array.
266:   */
267:   DMDAVecGetArray(da,X,&x);

269:   /*
270:      Get local grid boundaries (for 2-dimensional DMDA):
271:        xs, ys   - starting grid indices (no ghost points)
272:        xm, ym   - widths of local grid (no ghost points)

274:   */
275:   DMDAGetCorners(da,&xs,&ys,NULL,&xm,&ym,NULL);

277:   /*
278:      Compute initial guess over the locally owned part of the grid
279:   */
280:   for (j=ys; j<ys+ym; j++) {
281:     temp = (PetscReal)(PetscMin(j,My-j-1))*hy;
282:     for (i=xs; i<xs+xm; i++) {
283:       if (i == 0 || j == 0 || i == Mx-1 || j == My-1) {
284:         /* boundary conditions are all zero Dirichlet */
285:         x[j][i] = 0.0;
286:       } else {
287:         x[j][i] = temp1*PetscSqrtReal(PetscMin((PetscReal)(PetscMin(i,Mx-i-1))*hx,temp));
288:       }
289:     }
290:   }

292:   /*
293:      Restore vector
294:   */
295:   DMDAVecRestoreArray(da,X,&x);
296:   return(0);
297: }

299: /*
300:   FormExactSolution - Forms MMS solution

302:   Input Parameters:
303:   da - The DM
304:   user - user-defined application context

306:   Output Parameter:
307:   X - vector
308:  */
309: PetscErrorCode FormExactSolution(DM da, AppCtx *user, Vec U)
310: {
311:   DM             coordDA;
312:   Vec            coordinates;
313:   DMDACoor2d   **coords;
314:   PetscScalar  **u;
315:   PetscInt       xs, ys, xm, ym, i, j;

319:   DMDAGetCorners(da, &xs, &ys, NULL, &xm, &ym, NULL);
320:   DMGetCoordinateDM(da, &coordDA);
321:   DMGetCoordinates(da, &coordinates);
322:   DMDAVecGetArray(coordDA, coordinates, &coords);
323:   DMDAVecGetArray(da, U, &u);
324:   for (j = ys; j < ys+ym; ++j) {
325:     for (i = xs; i < xs+xm; ++i) {
326:       user->mms_solution(user,&coords[j][i],&u[j][i]);
327:     }
328:   }
329:   DMDAVecRestoreArray(da, U, &u);
330:   DMDAVecRestoreArray(coordDA, coordinates, &coords);
331:   return(0);
332: }

334: PetscErrorCode ZeroBCSolution(AppCtx *user,const DMDACoor2d *c,PetscScalar *u)
335: {
336:   u[0] = 0.;
337:   return 0;
338: }

340: /* The functions below evaluate the MMS solution u(x,y) and associated forcing

342:      f(x,y) = -u_xx - u_yy - lambda exp(u)

344:   such that u(x,y) is an exact solution with f(x,y) as the right hand side forcing term.
345:  */
346: PetscErrorCode MMSSolution1(AppCtx *user,const DMDACoor2d *c,PetscScalar *u)
347: {
348:   PetscReal x = PetscRealPart(c->x), y = PetscRealPart(c->y);
349:   u[0] = x*(1 - x)*y*(1 - y);
350:   PetscLogFlops(5);
351:   return 0;
352: }
353: PetscErrorCode MMSForcing1(AppCtx *user,const DMDACoor2d *c,PetscScalar *f)
354: {
355:   PetscReal x = PetscRealPart(c->x), y = PetscRealPart(c->y);
356:   f[0] = 2*x*(1 - x) + 2*y*(1 - y) - user->param*PetscExpReal(x*(1 - x)*y*(1 - y));
357:   return 0;
358: }

360: PetscErrorCode MMSSolution2(AppCtx *user,const DMDACoor2d *c,PetscScalar *u)
361: {
362:   PetscReal x = PetscRealPart(c->x), y = PetscRealPart(c->y);
363:   u[0] = PetscSinReal(PETSC_PI*x)*PetscSinReal(PETSC_PI*y);
364:   PetscLogFlops(5);
365:   return 0;
366: }
367: PetscErrorCode MMSForcing2(AppCtx *user,const DMDACoor2d *c,PetscScalar *f)
368: {
369:   PetscReal x = PetscRealPart(c->x), y = PetscRealPart(c->y);
370:   f[0] = 2*PetscSqr(PETSC_PI)*PetscSinReal(PETSC_PI*x)*PetscSinReal(PETSC_PI*y) - user->param*PetscExpReal(PetscSinReal(PETSC_PI*x)*PetscSinReal(PETSC_PI*y));
371:   return 0;
372: }

374: PetscErrorCode MMSSolution3(AppCtx *user,const DMDACoor2d *c,PetscScalar *u)
375: {
376:   PetscReal x = PetscRealPart(c->x), y = PetscRealPart(c->y);
377:   u[0] = PetscSinReal(user->m*PETSC_PI*x*(1-y))*PetscSinReal(user->n*PETSC_PI*y*(1-x));
378:   PetscLogFlops(5);
379:   return 0;
380: }
381: PetscErrorCode MMSForcing3(AppCtx *user,const DMDACoor2d *c,PetscScalar *f)
382: {
383:   PetscReal x = PetscRealPart(c->x), y = PetscRealPart(c->y);
384:   PetscReal m = user->m, n = user->n, lambda = user->param;
385:   f[0] = (-(PetscExpReal(PetscSinReal(m*PETSC_PI*x*(1 - y))*PetscSinReal(n*PETSC_PI*(1 - x)*y))*lambda)
386:           + PetscSqr(PETSC_PI)*(-2*m*n*((-1 + x)*x + (-1 + y)*y)*PetscCosReal(m*PETSC_PI*x*(-1 + y))*PetscCosReal(n*PETSC_PI*(-1 + x)*y)
387:                                 + (PetscSqr(m)*(PetscSqr(x) + PetscSqr(-1 + y)) + PetscSqr(n)*(PetscSqr(-1 + x) + PetscSqr(y)))
388:                                 *PetscSinReal(m*PETSC_PI*x*(-1 + y))*PetscSinReal(n*PETSC_PI*(-1 + x)*y)));
389:   return 0;
390: }

392: PetscErrorCode MMSSolution4(AppCtx *user,const DMDACoor2d *c,PetscScalar *u)
393: {
394:   const PetscReal Lx = 1.,Ly = 1.;
395:   PetscReal x = PetscRealPart(c->x), y = PetscRealPart(c->y);
396:   u[0] = (PetscPowReal(x,4)-PetscSqr(Lx)*PetscSqr(x))*(PetscPowReal(y,4)-PetscSqr(Ly)*PetscSqr(y));
397:   PetscLogFlops(9);
398:   return 0;
399: }
400: PetscErrorCode MMSForcing4(AppCtx *user,const DMDACoor2d *c,PetscScalar *f)
401: {
402:   const PetscReal Lx = 1.,Ly = 1.;
403:   PetscReal x = PetscRealPart(c->x), y = PetscRealPart(c->y);
404:   f[0] = (2*PetscSqr(x)*(PetscSqr(x)-PetscSqr(Lx))*(PetscSqr(Ly)-6*PetscSqr(y))
405:           + 2*PetscSqr(y)*(PetscSqr(Lx)-6*PetscSqr(x))*(PetscSqr(y)-PetscSqr(Ly))
406:           - user->param*PetscExpReal((PetscPowReal(x,4)-PetscSqr(Lx)*PetscSqr(x))*(PetscPowReal(y,4)-PetscSqr(Ly)*PetscSqr(y))));
407:   return 0;
408: }

410: /* ------------------------------------------------------------------- */
411: /*
412:    FormFunctionLocal - Evaluates nonlinear function, F(x) on local process patch

414:  */
415: PetscErrorCode FormFunctionLocal(DMDALocalInfo *info,PetscScalar **x,PetscScalar **f,AppCtx *user)
416: {
418:   PetscInt       i,j;
419:   PetscReal      lambda,hx,hy,hxdhy,hydhx;
420:   PetscScalar    u,ue,uw,un,us,uxx,uyy,mms_solution,mms_forcing;
421:   DMDACoor2d     c;

424:   lambda = user->param;
425:   hx     = 1.0/(PetscReal)(info->mx-1);
426:   hy     = 1.0/(PetscReal)(info->my-1);
427:   hxdhy  = hx/hy;
428:   hydhx  = hy/hx;
429:   /*
430:      Compute function over the locally owned part of the grid
431:   */
432:   for (j=info->ys; j<info->ys+info->ym; j++) {
433:     for (i=info->xs; i<info->xs+info->xm; i++) {
434:       if (i == 0 || j == 0 || i == info->mx-1 || j == info->my-1) {
435:         c.x = i*hx; c.y = j*hy;
436:         user->mms_solution(user,&c,&mms_solution);
437:         f[j][i] = 2.0*(hydhx+hxdhy)*(x[j][i] - mms_solution);
438:       } else {
439:         u  = x[j][i];
440:         uw = x[j][i-1];
441:         ue = x[j][i+1];
442:         un = x[j-1][i];
443:         us = x[j+1][i];

445:         /* Enforce boundary conditions at neighboring points -- setting these values causes the Jacobian to be symmetric. */
446:         if (i-1 == 0) {c.x = (i-1)*hx; c.y = j*hy; user->mms_solution(user,&c,&uw);}
447:         if (i+1 == info->mx-1) {c.x = (i+1)*hx; c.y = j*hy; user->mms_solution(user,&c,&ue);}
448:         if (j-1 == 0) {c.x = i*hx; c.y = (j-1)*hy; user->mms_solution(user,&c,&un);}
449:         if (j+1 == info->my-1) {c.x = i*hx; c.y = (j+1)*hy; user->mms_solution(user,&c,&us);}

451:         uxx     = (2.0*u - uw - ue)*hydhx;
452:         uyy     = (2.0*u - un - us)*hxdhy;
453:         mms_forcing = 0;
454:         c.x = i*hx; c.y = j*hy;
455:         if (user->mms_forcing) {user->mms_forcing(user,&c,&mms_forcing);}
456:         f[j][i] = uxx + uyy - hx*hy*(lambda*PetscExpScalar(u) + mms_forcing);
457:       }
458:     }
459:   }
460:   PetscLogFlops(11.0*info->ym*info->xm);
461:   return(0);
462: }

464: /* FormObjectiveLocal - Evaluates nonlinear function, F(x) on local process patch */
465: PetscErrorCode FormObjectiveLocal(DMDALocalInfo *info,PetscScalar **x,PetscReal *obj,AppCtx *user)
466: {
468:   PetscInt       i,j;
469:   PetscReal      lambda,hx,hy,hxdhy,hydhx,sc,lobj=0;
470:   PetscScalar    u,ue,uw,un,us,uxux,uyuy;
471:   MPI_Comm       comm;

474:   *obj   = 0;
475:   PetscObjectGetComm((PetscObject)info->da,&comm);
476:   lambda = user->param;
477:   hx     = 1.0/(PetscReal)(info->mx-1);
478:   hy     = 1.0/(PetscReal)(info->my-1);
479:   sc     = hx*hy*lambda;
480:   hxdhy  = hx/hy;
481:   hydhx  = hy/hx;
482:   /*
483:      Compute function over the locally owned part of the grid
484:   */
485:   for (j=info->ys; j<info->ys+info->ym; j++) {
486:     for (i=info->xs; i<info->xs+info->xm; i++) {
487:       if (i == 0 || j == 0 || i == info->mx-1 || j == info->my-1) {
488:         lobj += PetscRealPart((hydhx + hxdhy)*x[j][i]*x[j][i]);
489:       } else {
490:         u  = x[j][i];
491:         uw = x[j][i-1];
492:         ue = x[j][i+1];
493:         un = x[j-1][i];
494:         us = x[j+1][i];

496:         if (i-1 == 0) uw = 0.;
497:         if (i+1 == info->mx-1) ue = 0.;
498:         if (j-1 == 0) un = 0.;
499:         if (j+1 == info->my-1) us = 0.;

501:         /* F[u] = 1/2\int_{\omega}\nabla^2u(x)*u(x)*dx */

503:         uxux = u*(2.*u - ue - uw)*hydhx;
504:         uyuy = u*(2.*u - un - us)*hxdhy;

506:         lobj += PetscRealPart(0.5*(uxux + uyuy) - sc*PetscExpScalar(u));
507:       }
508:     }
509:   }
510:   PetscLogFlops(12.0*info->ym*info->xm);
511:   MPI_Allreduce(&lobj,obj,1,MPIU_REAL,MPIU_SUM,comm);
512:   return(0);
513: }

515: /*
516:    FormJacobianLocal - Evaluates Jacobian matrix on local process patch
517: */
518: PetscErrorCode FormJacobianLocal(DMDALocalInfo *info,PetscScalar **x,Mat jac,Mat jacpre,AppCtx *user)
519: {
521:   PetscInt       i,j,k;
522:   MatStencil     col[5],row;
523:   PetscScalar    lambda,v[5],hx,hy,hxdhy,hydhx,sc;
524:   DM             coordDA;
525:   Vec            coordinates;
526:   DMDACoor2d   **coords;

529:   lambda = user->param;
530:   /* Extract coordinates */
531:   DMGetCoordinateDM(info->da, &coordDA);
532:   DMGetCoordinates(info->da, &coordinates);
533:   DMDAVecGetArray(coordDA, coordinates, &coords);
534:   hx     = info->xm > 1 ? PetscRealPart(coords[info->ys][info->xs+1].x) - PetscRealPart(coords[info->ys][info->xs].x) : 1.0;
535:   hy     = info->ym > 1 ? PetscRealPart(coords[info->ys+1][info->xs].y) - PetscRealPart(coords[info->ys][info->xs].y) : 1.0;
536:   DMDAVecRestoreArray(coordDA, coordinates, &coords);
537:   hxdhy  = hx/hy;
538:   hydhx  = hy/hx;
539:   sc     = hx*hy*lambda;

541:   /*
542:      Compute entries for the locally owned part of the Jacobian.
543:       - Currently, all PETSc parallel matrix formats are partitioned by
544:         contiguous chunks of rows across the processors.
545:       - Each processor needs to insert only elements that it owns
546:         locally (but any non-local elements will be sent to the
547:         appropriate processor during matrix assembly).
548:       - Here, we set all entries for a particular row at once.
549:       - We can set matrix entries either using either
550:         MatSetValuesLocal() or MatSetValues(), as discussed above.
551:   */
552:   for (j=info->ys; j<info->ys+info->ym; j++) {
553:     for (i=info->xs; i<info->xs+info->xm; i++) {
554:       row.j = j; row.i = i;
555:       /* boundary points */
556:       if (i == 0 || j == 0 || i == info->mx-1 || j == info->my-1) {
557:         v[0] =  2.0*(hydhx + hxdhy);
558:         MatSetValuesStencil(jacpre,1,&row,1,&row,v,INSERT_VALUES);
559:       } else {
560:         k = 0;
561:         /* interior grid points */
562:         if (j-1 != 0) {
563:           v[k]     = -hxdhy;
564:           col[k].j = j - 1; col[k].i = i;
565:           k++;
566:         }
567:         if (i-1 != 0) {
568:           v[k]     = -hydhx;
569:           col[k].j = j;     col[k].i = i-1;
570:           k++;
571:         }

573:         v[k] = 2.0*(hydhx + hxdhy) - sc*PetscExpScalar(x[j][i]); col[k].j = row.j; col[k].i = row.i; k++;

575:         if (i+1 != info->mx-1) {
576:           v[k]     = -hydhx;
577:           col[k].j = j;     col[k].i = i+1;
578:           k++;
579:         }
580:         if (j+1 != info->mx-1) {
581:           v[k]     = -hxdhy;
582:           col[k].j = j + 1; col[k].i = i;
583:           k++;
584:         }
585:         MatSetValuesStencil(jacpre,1,&row,k,col,v,INSERT_VALUES);
586:       }
587:     }
588:   }

590:   /*
591:      Assemble matrix, using the 2-step process:
592:        MatAssemblyBegin(), MatAssemblyEnd().
593:   */
594:   MatAssemblyBegin(jacpre,MAT_FINAL_ASSEMBLY);
595:   MatAssemblyEnd(jacpre,MAT_FINAL_ASSEMBLY);
596:   /*
597:      Tell the matrix we will never add a new nonzero location to the
598:      matrix. If we do, it will generate an error.
599:   */
600:   MatSetOption(jac,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_TRUE);
601:   return(0);
602: }

604: PetscErrorCode FormFunctionMatlab(SNES snes,Vec X,Vec F,void *ptr)
605: {
606: #if PetscDefined(HAVE_MATLAB_ENGINE)
607:   AppCtx         *user = (AppCtx*)ptr;
609:   PetscInt       Mx,My;
610:   PetscReal      lambda,hx,hy;
611:   Vec            localX,localF;
612:   MPI_Comm       comm;
613:   DM             da;

616:   SNESGetDM(snes,&da);
617:   DMGetLocalVector(da,&localX);
618:   DMGetLocalVector(da,&localF);
619:   PetscObjectSetName((PetscObject)localX,"localX");
620:   PetscObjectSetName((PetscObject)localF,"localF");
621:   DMDAGetInfo(da,PETSC_IGNORE,&Mx,&My,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE);

623:   lambda = user->param;
624:   hx     = 1.0/(PetscReal)(Mx-1);
625:   hy     = 1.0/(PetscReal)(My-1);

627:   PetscObjectGetComm((PetscObject)snes,&comm);
628:   /*
629:      Scatter ghost points to local vector,using the 2-step process
630:         DMGlobalToLocalBegin(),DMGlobalToLocalEnd().
631:      By placing code between these two statements, computations can be
632:      done while messages are in transition.
633:   */
634:   DMGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
635:   DMGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
636:   PetscMatlabEnginePut(PETSC_MATLAB_ENGINE_(comm),(PetscObject)localX);
637:   PetscMatlabEngineEvaluate(PETSC_MATLAB_ENGINE_(comm),"localF=ex5m(localX,%18.16e,%18.16e,%18.16e)",hx,hy,lambda);
638:   PetscMatlabEngineGet(PETSC_MATLAB_ENGINE_(comm),(PetscObject)localF);

640:   /*
641:      Insert values into global vector
642:   */
643:   DMLocalToGlobalBegin(da,localF,INSERT_VALUES,F);
644:   DMLocalToGlobalEnd(da,localF,INSERT_VALUES,F);
645:   DMRestoreLocalVector(da,&localX);
646:   DMRestoreLocalVector(da,&localF);
647:   return(0);
648: #else
649:     return 0;                     /* Never called */
650: #endif
651: }

653: /* ------------------------------------------------------------------- */
654: /*
655:       Applies some sweeps on nonlinear Gauss-Seidel on each process

657:  */
658: PetscErrorCode NonlinearGS(SNES snes,Vec X, Vec B, void *ctx)
659: {
660:   PetscInt       i,j,k,Mx,My,xs,ys,xm,ym,its,tot_its,sweeps,l;
662:   PetscReal      lambda,hx,hy,hxdhy,hydhx,sc;
663:   PetscScalar    **x,**b,bij,F,F0=0,J,u,un,us,ue,eu,uw,uxx,uyy,y;
664:   PetscReal      atol,rtol,stol;
665:   DM             da;
666:   AppCtx         *user;
667:   Vec            localX,localB;

670:   tot_its = 0;
671:   SNESNGSGetSweeps(snes,&sweeps);
672:   SNESNGSGetTolerances(snes,&atol,&rtol,&stol,&its);
673:   SNESGetDM(snes,&da);
674:   DMGetApplicationContext(da,(void**)&user);

676:   DMDAGetInfo(da,PETSC_IGNORE,&Mx,&My,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE);

678:   lambda = user->param;
679:   hx     = 1.0/(PetscReal)(Mx-1);
680:   hy     = 1.0/(PetscReal)(My-1);
681:   sc     = hx*hy*lambda;
682:   hxdhy  = hx/hy;
683:   hydhx  = hy/hx;

685:   DMGetLocalVector(da,&localX);
686:   if (B) {
687:     DMGetLocalVector(da,&localB);
688:   }
689:   for (l=0; l<sweeps; l++) {

691:     DMGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
692:     DMGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
693:     if (B) {
694:       DMGlobalToLocalBegin(da,B,INSERT_VALUES,localB);
695:       DMGlobalToLocalEnd(da,B,INSERT_VALUES,localB);
696:     }
697:     /*
698:      Get a pointer to vector data.
699:      - For default PETSc vectors, VecGetArray() returns a pointer to
700:      the data array.  Otherwise, the routine is implementation dependent.
701:      - You MUST call VecRestoreArray() when you no longer need access to
702:      the array.
703:      */
704:     DMDAVecGetArray(da,localX,&x);
705:     if (B) DMDAVecGetArray(da,localB,&b);
706:     /*
707:      Get local grid boundaries (for 2-dimensional DMDA):
708:      xs, ys   - starting grid indices (no ghost points)
709:      xm, ym   - widths of local grid (no ghost points)
710:      */
711:     DMDAGetCorners(da,&xs,&ys,NULL,&xm,&ym,NULL);

713:     for (j=ys; j<ys+ym; j++) {
714:       for (i=xs; i<xs+xm; i++) {
715:         if (i == 0 || j == 0 || i == Mx-1 || j == My-1) {
716:           /* boundary conditions are all zero Dirichlet */
717:           x[j][i] = 0.0;
718:         } else {
719:           if (B) bij = b[j][i];
720:           else   bij = 0.;

722:           u  = x[j][i];
723:           un = x[j-1][i];
724:           us = x[j+1][i];
725:           ue = x[j][i-1];
726:           uw = x[j][i+1];

728:           for (k=0; k<its; k++) {
729:             eu  = PetscExpScalar(u);
730:             uxx = (2.0*u - ue - uw)*hydhx;
731:             uyy = (2.0*u - un - us)*hxdhy;
732:             F   = uxx + uyy - sc*eu - bij;
733:             if (k == 0) F0 = F;
734:             J  = 2.0*(hydhx + hxdhy) - sc*eu;
735:             y  = F/J;
736:             u -= y;
737:             tot_its++;

739:             if (atol > PetscAbsReal(PetscRealPart(F)) ||
740:                 rtol*PetscAbsReal(PetscRealPart(F0)) > PetscAbsReal(PetscRealPart(F)) ||
741:                 stol*PetscAbsReal(PetscRealPart(u)) > PetscAbsReal(PetscRealPart(y))) {
742:               break;
743:             }
744:           }
745:           x[j][i] = u;
746:         }
747:       }
748:     }
749:     /*
750:      Restore vector
751:      */
752:     DMDAVecRestoreArray(da,localX,&x);
753:     DMLocalToGlobalBegin(da,localX,INSERT_VALUES,X);
754:     DMLocalToGlobalEnd(da,localX,INSERT_VALUES,X);
755:   }
756:   PetscLogFlops(tot_its*(21.0));
757:   DMRestoreLocalVector(da,&localX);
758:   if (B) {
759:     DMDAVecRestoreArray(da,localB,&b);
760:     DMRestoreLocalVector(da,&localB);
761:   }
762:   return(0);
763: }

765: /*TEST

767:    test:
768:      suffix: asm_0
769:      requires: !single
770:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 2 -pc_asm_overlap 0 -pc_asm_local_type additive -sub_pc_type lu

772:    test:
773:      suffix: msm_0
774:      requires: !single
775:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 2 -pc_asm_overlap 0 -pc_asm_local_type multiplicative -sub_pc_type lu

777:    test:
778:      suffix: asm_1
779:      requires: !single
780:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 2 -pc_asm_overlap 0 -pc_asm_local_type additive -sub_pc_type lu -da_grid_x 8

782:    test:
783:      suffix: msm_1
784:      requires: !single
785:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 2 -pc_asm_overlap 0 -pc_asm_local_type multiplicative -sub_pc_type lu -da_grid_x 8

787:    test:
788:      suffix: asm_2
789:      requires: !single
790:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 3 -pc_asm_overlap 0 -pc_asm_local_type additive -sub_pc_type lu -da_grid_x 8

792:    test:
793:      suffix: msm_2
794:      requires: !single
795:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 3 -pc_asm_overlap 0 -pc_asm_local_type multiplicative -sub_pc_type lu -da_grid_x 8

797:    test:
798:      suffix: asm_3
799:      requires: !single
800:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 4 -pc_asm_overlap 0 -pc_asm_local_type additive -sub_pc_type lu -da_grid_x 8

802:    test:
803:      suffix: msm_3
804:      requires: !single
805:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 4 -pc_asm_overlap 0 -pc_asm_local_type multiplicative -sub_pc_type lu -da_grid_x 8

807:    test:
808:      suffix: asm_4
809:      requires: !single
810:      nsize: 2
811:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 2 -pc_asm_overlap 0 -pc_asm_local_type additive -sub_pc_type lu -da_grid_x 8

813:    test:
814:      suffix: msm_4
815:      requires: !single
816:      nsize: 2
817:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 2 -pc_asm_overlap 0 -pc_asm_local_type multiplicative -sub_pc_type lu -da_grid_x 8

819:    test:
820:      suffix: asm_5
821:      requires: !single
822:      nsize: 2
823:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 4 -pc_asm_overlap 0 -pc_asm_local_type additive -sub_pc_type lu -da_grid_x 8

825:    test:
826:      suffix: msm_5
827:      requires: !single
828:      nsize: 2
829:      args: -mms 1 -par 0.0 -snes_monitor_short -snes_converged_reason -snes_view -ksp_rtol 1.0e-9 -ksp_monitor_short -ksp_type richardson -pc_type asm -pc_asm_blocks 4 -pc_asm_overlap 0 -pc_asm_local_type multiplicative -sub_pc_type lu -da_grid_x 8

831:    test:
832:      args: -snes_rtol 1.e-5 -pc_type mg -ksp_monitor_short -snes_view -pc_mg_levels 3 -pc_mg_galerkin pmat -da_grid_x 17 -da_grid_y 17 -mg_levels_ksp_monitor_short -mg_levels_ksp_norm_type unpreconditioned -snes_monitor_short -mg_levels_ksp_chebyshev_esteig 0.5,1.1 -mg_levels_pc_type sor -pc_mg_type full
833:      requires: !single

835:    test:
836:      suffix: 2
837:      args: -pc_type mg -ksp_converged_reason -snes_view -pc_mg_galerkin pmat -snes_grid_sequence 3 -mg_levels_ksp_norm_type unpreconditioned -snes_monitor_short -mg_levels_ksp_chebyshev_esteig 0.5,1.1 -mg_levels_pc_type sor -pc_mg_type full -ksp_atol -1.

839:    test:
840:      suffix: 3
841:      nsize: 2
842:      args: -snes_grid_sequence 2 -snes_mf_operator -snes_converged_reason -snes_view -pc_type mg -snes_atol -1 -snes_rtol 1.e-2
843:      filter: grep -v "otal number of function evaluations"

845:    test:
846:      suffix: 4
847:      nsize: 2
848:      args: -snes_grid_sequence 2 -snes_monitor_short -ksp_converged_reason -snes_converged_reason -snes_view -pc_type mg -snes_atol -1 -ksp_atol -1

850:    test:
851:      suffix: 5_anderson
852:      args: -da_grid_x 81 -da_grid_y 81 -snes_monitor_short -snes_max_it 50 -par 6.0 -snes_type anderson

854:    test:
855:      suffix: 5_aspin
856:      nsize: 4
857:      args: -snes_monitor_short -ksp_monitor_short -snes_converged_reason -da_refine 4 -da_overlap 3 -snes_type aspin -snes_view

859:    test:
860:      suffix: 5_broyden
861:      args: -da_grid_x 81 -da_grid_y 81 -snes_monitor_short -snes_max_it 50 -par 6.0 -snes_type qn -snes_qn_type broyden -snes_qn_m 10

863:    test:
864:      suffix: 5_fas
865:      args: -fas_coarse_snes_max_it 1 -fas_coarse_pc_type lu -fas_coarse_ksp_type preonly -snes_monitor_short -snes_type fas -fas_coarse_ksp_type richardson -da_refine 6
866:      requires: !single

868:    test:
869:      suffix: 5_fas_additive
870:      args: -fas_coarse_snes_max_it 1 -fas_coarse_pc_type lu -fas_coarse_ksp_type preonly -snes_monitor_short -snes_type fas -fas_coarse_ksp_type richardson -da_refine 6 -snes_fas_type additive -snes_max_it 50

872:    test:
873:      suffix: 5_fas_monitor
874:      args: -da_refine 1 -snes_type fas -snes_fas_monitor
875:      requires: !single

877:    test:
878:      suffix: 5_ls
879:      args: -da_grid_x 81 -da_grid_y 81 -snes_monitor_short -snes_max_it 50 -par 6.0 -snes_type newtonls

881:    test:
882:      suffix: 5_ls_sell_sor
883:      args: -da_grid_x 81 -da_grid_y 81 -snes_monitor_short -snes_max_it 50 -par 6.0 -snes_type newtonls -dm_mat_type sell -pc_type sor
884:      output_file: output/ex5_5_ls.out

886:    test:
887:      suffix: 5_nasm
888:      nsize: 4
889:      args: -snes_monitor_short -snes_converged_reason -da_refine 4 -da_overlap 3 -snes_type nasm -snes_nasm_type restrict -snes_max_it 10

891:    test:
892:      suffix: 5_ncg
893:      args: -da_grid_x 81 -da_grid_y 81 -snes_monitor_short -snes_max_it 50 -par 6.0 -snes_type ncg -snes_ncg_type fr

895:    test:
896:      suffix: 5_newton_asm_dmda
897:      nsize: 4
898:      args: -snes_monitor_short -ksp_monitor_short -snes_converged_reason -da_refine 4 -da_overlap 3 -snes_type newtonls -pc_type asm -pc_asm_dm_subdomains -malloc_dump
899:      requires: !single

901:    test:
902:      suffix: 5_newton_gasm_dmda
903:      nsize: 4
904:      args: -snes_monitor_short -ksp_monitor_short -snes_converged_reason -da_refine 4 -da_overlap 3 -snes_type newtonls -pc_type gasm -malloc_dump
905:      requires: !single

907:    test:
908:      suffix: 5_ngmres
909:      args: -da_grid_x 81 -da_grid_y 81 -snes_monitor_short -snes_max_it 50 -par 6.0 -snes_type ngmres -snes_ngmres_m 10

911:    test:
912:      suffix: 5_ngmres_fas
913:      args: -snes_rtol 1.e-4 -snes_type ngmres -npc_fas_coarse_snes_max_it 1 -npc_fas_coarse_snes_type newtonls -npc_fas_coarse_pc_type lu -npc_fas_coarse_ksp_type preonly -snes_ngmres_m 10 -snes_monitor_short -npc_snes_max_it 1 -npc_snes_type fas -npc_fas_coarse_ksp_type richardson -da_refine 6

915:    test:
916:      suffix: 5_ngmres_ngs
917:      args: -da_grid_x 81 -da_grid_y 81 -snes_monitor_short -snes_max_it 50 -par 6.0 -snes_type ngmres -npc_snes_type ngs -npc_snes_max_it 1

919:    test:
920:      suffix: 5_ngmres_nrichardson
921:      args: -da_grid_x 81 -da_grid_y 81 -snes_monitor_short -snes_max_it 50 -par 6.0 -snes_type ngmres -snes_ngmres_m 10 -npc_snes_type nrichardson -npc_snes_max_it 3
922:      output_file: output/ex5_5_ngmres_richardson.out

924:    test:
925:      suffix: 5_nrichardson
926:      args: -da_grid_x 81 -da_grid_y 81 -snes_monitor_short -snes_max_it 50 -par 6.0 -snes_type nrichardson

928:    test:
929:      suffix: 5_qn
930:      args: -da_grid_x 81 -da_grid_y 81 -snes_monitor_short -snes_max_it 50 -par 6.0 -snes_type qn -snes_linesearch_type cp -snes_qn_m 10

932:    test:
933:      suffix: 6
934:      nsize: 4
935:      args: -snes_converged_reason -ksp_converged_reason -da_grid_x 129 -da_grid_y 129 -pc_type mg -pc_mg_levels 8 -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_esteig 0,0.5,0,1.1 -mg_levels_ksp_max_it 2

937:    test:
938:      requires: complex !single
939:      suffix: complex
940:      args: -snes_mf_operator -mat_mffd_complex -snes_monitor

942: TEST*/