Actual source code: ex3.c
petsc-3.9.4 2018-09-11
1: static char help[] = "Newton methods to solve u'' + u^{2} = f in parallel.\n\
2: This example employs a user-defined monitoring routine and optionally a user-defined\n\
3: routine to check candidate iterates produced by line search routines.\n\
4: The command line options include:\n\
5: -pre_check_iterates : activate checking of iterates\n\
6: -post_check_iterates : activate checking of iterates\n\
7: -check_tol <tol>: set tolerance for iterate checking\n\
8: -user_precond : activate a (trivial) user-defined preconditioner\n\n";
10: /*T
11: Concepts: SNES^basic parallel example
12: Concepts: SNES^setting a user-defined monitoring routine
13: Processors: n
14: T*/
18: /*
19: Include "petscdraw.h" so that we can use distributed arrays (DMDAs).
20: Include "petscdraw.h" so that we can use PETSc drawing routines.
21: Include "petscsnes.h" so that we can use SNES solvers. Note that this
22: file automatically includes:
23: petscsys.h - base PETSc routines
24: petscvec.h - vectors
25: petscmat.h - matrices
26: petscis.h - index sets
27: petscksp.h - Krylov subspace methods
28: petscviewer.h - viewers
29: petscpc.h - preconditioners
30: petscksp.h - linear solvers
31: */
33: #include <petscdm.h>
34: #include <petscdmda.h>
35: #include <petscsnes.h>
37: /*
38: User-defined routines.
39: */
40: PetscErrorCode FormJacobian(SNES,Vec,Mat,Mat,void*);
41: PetscErrorCode FormFunction(SNES,Vec,Vec,void*);
42: PetscErrorCode FormInitialGuess(Vec);
43: PetscErrorCode Monitor(SNES,PetscInt,PetscReal,void*);
44: PetscErrorCode PreCheck(SNESLineSearch,Vec,Vec,PetscBool*,void*);
45: PetscErrorCode PostCheck(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*);
46: PetscErrorCode PostSetSubKSP(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*);
47: PetscErrorCode MatrixFreePreconditioner(PC,Vec,Vec);
49: /*
50: User-defined application context
51: */
52: typedef struct {
53: DM da; /* distributed array */
54: Vec F; /* right-hand-side of PDE */
55: PetscMPIInt rank; /* rank of processor */
56: PetscMPIInt size; /* size of communicator */
57: PetscReal h; /* mesh spacing */
58: } ApplicationCtx;
60: /*
61: User-defined context for monitoring
62: */
63: typedef struct {
64: PetscViewer viewer;
65: } MonitorCtx;
67: /*
68: User-defined context for checking candidate iterates that are
69: determined by line search methods
70: */
71: typedef struct {
72: Vec last_step; /* previous iterate */
73: PetscReal tolerance; /* tolerance for changes between successive iterates */
74: ApplicationCtx *user;
75: } StepCheckCtx;
77: typedef struct {
78: PetscInt its0; /* num of prevous outer KSP iterations */
79: } SetSubKSPCtx;
81: int main(int argc,char **argv)
82: {
83: SNES snes; /* SNES context */
84: SNESLineSearch linesearch; /* SNESLineSearch context */
85: Mat J; /* Jacobian matrix */
86: ApplicationCtx ctx; /* user-defined context */
87: Vec x,r,U,F; /* vectors */
88: MonitorCtx monP; /* monitoring context */
89: StepCheckCtx checkP; /* step-checking context */
90: SetSubKSPCtx checkP1;
91: PetscBool pre_check,post_check,post_setsubksp; /* flag indicating whether we're checking candidate iterates */
92: PetscScalar xp,*FF,*UU,none = -1.0;
94: PetscInt its,N = 5,i,maxit,maxf,xs,xm;
95: PetscReal abstol,rtol,stol,norm;
96: PetscBool flg;
98: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
99: MPI_Comm_rank(PETSC_COMM_WORLD,&ctx.rank);
100: MPI_Comm_size(PETSC_COMM_WORLD,&ctx.size);
101: PetscOptionsGetInt(NULL,NULL,"-n",&N,NULL);
102: ctx.h = 1.0/(N-1);
104: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
105: Create nonlinear solver context
106: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
108: SNESCreate(PETSC_COMM_WORLD,&snes);
110: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
111: Create vector data structures; set function evaluation routine
112: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
114: /*
115: Create distributed array (DMDA) to manage parallel grid and vectors
116: */
117: DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,N,1,1,NULL,&ctx.da);
118: DMSetFromOptions(ctx.da);
119: DMSetUp(ctx.da);
121: /*
122: Extract global and local vectors from DMDA; then duplicate for remaining
123: vectors that are the same types
124: */
125: DMCreateGlobalVector(ctx.da,&x);
126: VecDuplicate(x,&r);
127: VecDuplicate(x,&F); ctx.F = F;
128: VecDuplicate(x,&U);
130: /*
131: Set function evaluation routine and vector. Whenever the nonlinear
132: solver needs to compute the nonlinear function, it will call this
133: routine.
134: - Note that the final routine argument is the user-defined
135: context that provides application-specific data for the
136: function evaluation routine.
137: */
138: SNESSetFunction(snes,r,FormFunction,&ctx);
140: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
141: Create matrix data structure; set Jacobian evaluation routine
142: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
144: MatCreate(PETSC_COMM_WORLD,&J);
145: MatSetSizes(J,PETSC_DECIDE,PETSC_DECIDE,N,N);
146: MatSetFromOptions(J);
147: MatSeqAIJSetPreallocation(J,3,NULL);
148: MatMPIAIJSetPreallocation(J,3,NULL,3,NULL);
150: /*
151: Set Jacobian matrix data structure and default Jacobian evaluation
152: routine. Whenever the nonlinear solver needs to compute the
153: Jacobian matrix, it will call this routine.
154: - Note that the final routine argument is the user-defined
155: context that provides application-specific data for the
156: Jacobian evaluation routine.
157: */
158: SNESSetJacobian(snes,J,J,FormJacobian,&ctx);
160: /*
161: Optionally allow user-provided preconditioner
162: */
163: PetscOptionsHasName(NULL,NULL,"-user_precond",&flg);
164: if (flg) {
165: KSP ksp;
166: PC pc;
167: SNESGetKSP(snes,&ksp);
168: KSPGetPC(ksp,&pc);
169: PCSetType(pc,PCSHELL);
170: PCShellSetApply(pc,MatrixFreePreconditioner);
171: }
173: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
174: Customize nonlinear solver; set runtime options
175: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
177: /*
178: Set an optional user-defined monitoring routine
179: */
180: PetscViewerDrawOpen(PETSC_COMM_WORLD,0,0,0,0,400,400,&monP.viewer);
181: SNESMonitorSet(snes,Monitor,&monP,0);
183: /*
184: Set names for some vectors to facilitate monitoring (optional)
185: */
186: PetscObjectSetName((PetscObject)x,"Approximate Solution");
187: PetscObjectSetName((PetscObject)U,"Exact Solution");
189: /*
190: Set SNES/KSP/KSP/PC runtime options, e.g.,
191: -snes_view -snes_monitor -ksp_type <ksp> -pc_type <pc>
192: */
193: SNESSetFromOptions(snes);
195: /*
196: Set an optional user-defined routine to check the validity of candidate
197: iterates that are determined by line search methods
198: */
199: SNESGetLineSearch(snes, &linesearch);
200: PetscOptionsHasName(NULL,NULL,"-post_check_iterates",&post_check);
202: if (post_check) {
203: PetscPrintf(PETSC_COMM_WORLD,"Activating post step checking routine\n");
204: SNESLineSearchSetPostCheck(linesearch,PostCheck,&checkP);
205: VecDuplicate(x,&(checkP.last_step));
207: checkP.tolerance = 1.0;
208: checkP.user = &ctx;
210: PetscOptionsGetReal(NULL,NULL,"-check_tol",&checkP.tolerance,NULL);
211: }
213: PetscOptionsHasName(NULL,NULL,"-post_setsubksp",&post_setsubksp);
214: if (post_setsubksp) {
215: PetscPrintf(PETSC_COMM_WORLD,"Activating post setsubksp\n");
216: SNESLineSearchSetPostCheck(linesearch,PostSetSubKSP,&checkP1);
217: }
219: PetscOptionsHasName(NULL,NULL,"-pre_check_iterates",&pre_check);
220: if (pre_check) {
221: PetscPrintf(PETSC_COMM_WORLD,"Activating pre step checking routine\n");
222: SNESLineSearchSetPreCheck(linesearch,PreCheck,&checkP);
223: }
225: /*
226: Print parameters used for convergence testing (optional) ... just
227: to demonstrate this routine; this information is also printed with
228: the option -snes_view
229: */
230: SNESGetTolerances(snes,&abstol,&rtol,&stol,&maxit,&maxf);
231: PetscPrintf(PETSC_COMM_WORLD,"atol=%g, rtol=%g, stol=%g, maxit=%D, maxf=%D\n",(double)abstol,(double)rtol,(double)stol,maxit,maxf);
233: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
234: Initialize application:
235: Store right-hand-side of PDE and exact solution
236: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
238: /*
239: Get local grid boundaries (for 1-dimensional DMDA):
240: xs, xm - starting grid index, width of local grid (no ghost points)
241: */
242: DMDAGetCorners(ctx.da,&xs,NULL,NULL,&xm,NULL,NULL);
244: /*
245: Get pointers to vector data
246: */
247: DMDAVecGetArray(ctx.da,F,&FF);
248: DMDAVecGetArray(ctx.da,U,&UU);
250: /*
251: Compute local vector entries
252: */
253: xp = ctx.h*xs;
254: for (i=xs; i<xs+xm; i++) {
255: FF[i] = 6.0*xp + PetscPowScalar(xp+1.e-12,6.0); /* +1.e-12 is to prevent 0^6 */
256: UU[i] = xp*xp*xp;
257: xp += ctx.h;
258: }
260: /*
261: Restore vectors
262: */
263: DMDAVecRestoreArray(ctx.da,F,&FF);
264: DMDAVecRestoreArray(ctx.da,U,&UU);
266: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
267: Evaluate initial guess; then solve nonlinear system
268: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
270: /*
271: Note: The user should initialize the vector, x, with the initial guess
272: for the nonlinear solver prior to calling SNESSolve(). In particular,
273: to employ an initial guess of zero, the user should explicitly set
274: this vector to zero by calling VecSet().
275: */
276: FormInitialGuess(x);
277: SNESSolve(snes,NULL,x);
278: SNESGetIterationNumber(snes,&its);
279: PetscPrintf(PETSC_COMM_WORLD,"Number of SNES iterations = %D\n",its);
281: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
282: Check solution and clean up
283: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
285: /*
286: Check the error
287: */
288: VecAXPY(x,none,U);
289: VecNorm(x,NORM_2,&norm);
290: PetscPrintf(PETSC_COMM_WORLD,"Norm of error %g Iterations %D\n",(double)norm,its);
292: /*
293: Free work space. All PETSc objects should be destroyed when they
294: are no longer needed.
295: */
296: PetscViewerDestroy(&monP.viewer);
297: if (post_check) {VecDestroy(&checkP.last_step);}
298: VecDestroy(&x);
299: VecDestroy(&r);
300: VecDestroy(&U);
301: VecDestroy(&F);
302: MatDestroy(&J);
303: SNESDestroy(&snes);
304: DMDestroy(&ctx.da);
305: PetscFinalize();
306: return ierr;
307: }
309: /* ------------------------------------------------------------------- */
310: /*
311: FormInitialGuess - Computes initial guess.
313: Input/Output Parameter:
314: . x - the solution vector
315: */
316: PetscErrorCode FormInitialGuess(Vec x)
317: {
319: PetscScalar pfive = .50;
322: VecSet(x,pfive);
323: return(0);
324: }
326: /* ------------------------------------------------------------------- */
327: /*
328: FormFunction - Evaluates nonlinear function, F(x).
330: Input Parameters:
331: . snes - the SNES context
332: . x - input vector
333: . ctx - optional user-defined context, as set by SNESSetFunction()
335: Output Parameter:
336: . f - function vector
338: Note:
339: The user-defined context can contain any application-specific
340: data needed for the function evaluation.
341: */
342: PetscErrorCode FormFunction(SNES snes,Vec x,Vec f,void *ctx)
343: {
344: ApplicationCtx *user = (ApplicationCtx*) ctx;
345: DM da = user->da;
346: PetscScalar *xx,*ff,*FF,d;
348: PetscInt i,M,xs,xm;
349: Vec xlocal;
352: DMGetLocalVector(da,&xlocal);
353: /*
354: Scatter ghost points to local vector, using the 2-step process
355: DMGlobalToLocalBegin(), DMGlobalToLocalEnd().
356: By placing code between these two statements, computations can
357: be done while messages are in transition.
358: */
359: DMGlobalToLocalBegin(da,x,INSERT_VALUES,xlocal);
360: DMGlobalToLocalEnd(da,x,INSERT_VALUES,xlocal);
362: /*
363: Get pointers to vector data.
364: - The vector xlocal includes ghost point; the vectors x and f do
365: NOT include ghost points.
366: - Using DMDAVecGetArray() allows accessing the values using global ordering
367: */
368: DMDAVecGetArray(da,xlocal,&xx);
369: DMDAVecGetArray(da,f,&ff);
370: DMDAVecGetArray(da,user->F,&FF);
372: /*
373: Get local grid boundaries (for 1-dimensional DMDA):
374: xs, xm - starting grid index, width of local grid (no ghost points)
375: */
376: DMDAGetCorners(da,&xs,NULL,NULL,&xm,NULL,NULL);
377: DMDAGetInfo(da,NULL,&M,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);
379: /*
380: Set function values for boundary points; define local interior grid point range:
381: xsi - starting interior grid index
382: xei - ending interior grid index
383: */
384: if (xs == 0) { /* left boundary */
385: ff[0] = xx[0];
386: xs++;xm--;
387: }
388: if (xs+xm == M) { /* right boundary */
389: ff[xs+xm-1] = xx[xs+xm-1] - 1.0;
390: xm--;
391: }
393: /*
394: Compute function over locally owned part of the grid (interior points only)
395: */
396: d = 1.0/(user->h*user->h);
397: for (i=xs; i<xs+xm; i++) ff[i] = d*(xx[i-1] - 2.0*xx[i] + xx[i+1]) + xx[i]*xx[i] - FF[i];
399: /*
400: Restore vectors
401: */
402: DMDAVecRestoreArray(da,xlocal,&xx);
403: DMDAVecRestoreArray(da,f,&ff);
404: DMDAVecRestoreArray(da,user->F,&FF);
405: DMRestoreLocalVector(da,&xlocal);
406: return(0);
407: }
409: /* ------------------------------------------------------------------- */
410: /*
411: FormJacobian - Evaluates Jacobian matrix.
413: Input Parameters:
414: . snes - the SNES context
415: . x - input vector
416: . dummy - optional user-defined context (not used here)
418: Output Parameters:
419: . jac - Jacobian matrix
420: . B - optionally different preconditioning matrix
421: . flag - flag indicating matrix structure
422: */
423: PetscErrorCode FormJacobian(SNES snes,Vec x,Mat jac,Mat B,void *ctx)
424: {
425: ApplicationCtx *user = (ApplicationCtx*) ctx;
426: PetscScalar *xx,d,A[3];
428: PetscInt i,j[3],M,xs,xm;
429: DM da = user->da;
432: /*
433: Get pointer to vector data
434: */
435: DMDAVecGetArrayRead(da,x,&xx);
436: DMDAGetCorners(da,&xs,NULL,NULL,&xm,NULL,NULL);
438: /*
439: Get range of locally owned matrix
440: */
441: DMDAGetInfo(da,NULL,&M,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);
443: /*
444: Determine starting and ending local indices for interior grid points.
445: Set Jacobian entries for boundary points.
446: */
448: if (xs == 0) { /* left boundary */
449: i = 0; A[0] = 1.0;
451: MatSetValues(jac,1,&i,1,&i,A,INSERT_VALUES);
452: xs++;xm--;
453: }
454: if (xs+xm == M) { /* right boundary */
455: i = M-1;
456: A[0] = 1.0;
457: MatSetValues(jac,1,&i,1,&i,A,INSERT_VALUES);
458: xm--;
459: }
461: /*
462: Interior grid points
463: - Note that in this case we set all elements for a particular
464: row at once.
465: */
466: d = 1.0/(user->h*user->h);
467: for (i=xs; i<xs+xm; i++) {
468: j[0] = i - 1; j[1] = i; j[2] = i + 1;
469: A[0] = A[2] = d; A[1] = -2.0*d + 2.0*xx[i];
470: MatSetValues(jac,1,&i,3,j,A,INSERT_VALUES);
471: }
473: /*
474: Assemble matrix, using the 2-step process:
475: MatAssemblyBegin(), MatAssemblyEnd().
476: By placing code between these two statements, computations can be
477: done while messages are in transition.
479: Also, restore vector.
480: */
482: MatAssemblyBegin(jac,MAT_FINAL_ASSEMBLY);
483: DMDAVecRestoreArrayRead(da,x,&xx);
484: MatAssemblyEnd(jac,MAT_FINAL_ASSEMBLY);
486: return(0);
487: }
489: /* ------------------------------------------------------------------- */
490: /*
491: Monitor - Optional user-defined monitoring routine that views the
492: current iterate with an x-window plot. Set by SNESMonitorSet().
494: Input Parameters:
495: snes - the SNES context
496: its - iteration number
497: norm - 2-norm function value (may be estimated)
498: ctx - optional user-defined context for private data for the
499: monitor routine, as set by SNESMonitorSet()
501: Note:
502: See the manpage for PetscViewerDrawOpen() for useful runtime options,
503: such as -nox to deactivate all x-window output.
504: */
505: PetscErrorCode Monitor(SNES snes,PetscInt its,PetscReal fnorm,void *ctx)
506: {
508: MonitorCtx *monP = (MonitorCtx*) ctx;
509: Vec x;
512: PetscPrintf(PETSC_COMM_WORLD,"iter = %D,SNES Function norm %g\n",its,(double)fnorm);
513: SNESGetSolution(snes,&x);
514: VecView(x,monP->viewer);
515: return(0);
516: }
518: /* ------------------------------------------------------------------- */
519: /*
520: PreCheck - Optional user-defined routine that checks the validity of
521: candidate steps of a line search method. Set by SNESLineSearchSetPreCheck().
523: Input Parameters:
524: snes - the SNES context
525: xcurrent - current solution
526: y - search direction and length
528: Output Parameters:
529: y - proposed step (search direction and length) (possibly changed)
530: changed_y - tells if the step has changed or not
531: */
532: PetscErrorCode PreCheck(SNESLineSearch linesearch,Vec xcurrent,Vec y, PetscBool *changed_y, void * ctx)
533: {
535: *changed_y = PETSC_FALSE;
536: return(0);
537: }
539: /* ------------------------------------------------------------------- */
540: /*
541: PostCheck - Optional user-defined routine that checks the validity of
542: candidate steps of a line search method. Set by SNESLineSearchSetPostCheck().
544: Input Parameters:
545: snes - the SNES context
546: ctx - optional user-defined context for private data for the
547: monitor routine, as set by SNESLineSearchSetPostCheck()
548: xcurrent - current solution
549: y - search direction and length
550: x - the new candidate iterate
552: Output Parameters:
553: y - proposed step (search direction and length) (possibly changed)
554: x - current iterate (possibly modified)
556: */
557: PetscErrorCode PostCheck(SNESLineSearch linesearch,Vec xcurrent,Vec y,Vec x,PetscBool *changed_y,PetscBool *changed_x, void * ctx)
558: {
560: PetscInt i,iter,xs,xm;
561: StepCheckCtx *check;
562: ApplicationCtx *user;
563: PetscScalar *xa,*xa_last,tmp;
564: PetscReal rdiff;
565: DM da;
566: SNES snes;
569: *changed_x = PETSC_FALSE;
570: *changed_y = PETSC_FALSE;
572: SNESLineSearchGetSNES(linesearch, &snes);
573: check = (StepCheckCtx*)ctx;
574: user = check->user;
575: SNESGetIterationNumber(snes,&iter);
577: /* iteration 1 indicates we are working on the second iteration */
578: if (iter > 0) {
579: da = user->da;
580: PetscPrintf(PETSC_COMM_WORLD,"Checking candidate step at iteration %D with tolerance %g\n",iter,(double)check->tolerance);
582: /* Access local array data */
583: DMDAVecGetArray(da,check->last_step,&xa_last);
584: DMDAVecGetArray(da,x,&xa);
585: DMDAGetCorners(da,&xs,NULL,NULL,&xm,NULL,NULL);
587: /*
588: If we fail the user-defined check for validity of the candidate iterate,
589: then modify the iterate as we like. (Note that the particular modification
590: below is intended simply to demonstrate how to manipulate this data, not
591: as a meaningful or appropriate choice.)
592: */
593: for (i=xs; i<xs+xm; i++) {
594: if (!PetscAbsScalar(xa[i])) rdiff = 2*check->tolerance;
595: else rdiff = PetscAbsScalar((xa[i] - xa_last[i])/xa[i]);
596: if (rdiff > check->tolerance) {
597: tmp = xa[i];
598: xa[i] = .5*(xa[i] + xa_last[i]);
599: *changed_x = PETSC_TRUE;
600: PetscPrintf(PETSC_COMM_WORLD," Altering entry %D: x=%g, x_last=%g, diff=%g, x_new=%g\n",i,(double)PetscAbsScalar(tmp),(double)PetscAbsScalar(xa_last[i]),(double)rdiff,(double)PetscAbsScalar(xa[i]));
601: }
602: }
603: DMDAVecRestoreArray(da,check->last_step,&xa_last);
604: DMDAVecRestoreArray(da,x,&xa);
605: }
606: VecCopy(x,check->last_step);
607: return(0);
608: }
610: /* ------------------------------------------------------------------- */
611: /*
612: PostSetSubKSP - Optional user-defined routine that reset SubKSP options when hierarchical bjacobi PC is used
613: e.g,
614: mpiexec -n 8 ./ex3 -nox -n 10000 -ksp_type fgmres -pc_type bjacobi -pc_bjacobi_blocks 4 -sub_ksp_type gmres -sub_ksp_max_it 3 -post_setsubksp -sub_ksp_rtol 1.e-16
615: Set by SNESLineSearchSetPostCheck().
617: Input Parameters:
618: linesearch - the LineSearch context
619: xcurrent - current solution
620: y - search direction and length
621: x - the new candidate iterate
623: Output Parameters:
624: y - proposed step (search direction and length) (possibly changed)
625: x - current iterate (possibly modified)
627: */
628: PetscErrorCode PostSetSubKSP(SNESLineSearch linesearch,Vec xcurrent,Vec y,Vec x,PetscBool *changed_y,PetscBool *changed_x, void * ctx)
629: {
631: SetSubKSPCtx *check;
632: PetscInt iter,its,sub_its,maxit;
633: KSP ksp,sub_ksp,*sub_ksps;
634: PC pc;
635: PetscReal ksp_ratio;
636: SNES snes;
639: SNESLineSearchGetSNES(linesearch, &snes);
640: check = (SetSubKSPCtx*)ctx;
641: SNESGetIterationNumber(snes,&iter);
642: SNESGetKSP(snes,&ksp);
643: KSPGetPC(ksp,&pc);
644: PCBJacobiGetSubKSP(pc,NULL,NULL,&sub_ksps);
645: sub_ksp = sub_ksps[0];
646: KSPGetIterationNumber(ksp,&its); /* outer KSP iteration number */
647: KSPGetIterationNumber(sub_ksp,&sub_its); /* inner KSP iteration number */
649: if (iter) {
650: PetscPrintf(PETSC_COMM_WORLD," ...PostCheck snes iteration %D, ksp_it %d %d, subksp_it %d\n",iter,check->its0,its,sub_its);
651: ksp_ratio = ((PetscReal)(its))/check->its0;
652: maxit = (PetscInt)(ksp_ratio*sub_its + 0.5);
653: if (maxit < 2) maxit = 2;
654: KSPSetTolerances(sub_ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,maxit);
655: PetscPrintf(PETSC_COMM_WORLD," ...ksp_ratio %g, new maxit %d\n\n",ksp_ratio,maxit);
656: }
657: check->its0 = its; /* save current outer KSP iteration number */
658: return(0);
659: }
661: /* ------------------------------------------------------------------- */
662: /*
663: MatrixFreePreconditioner - This routine demonstrates the use of a
664: user-provided preconditioner. This code implements just the null
665: preconditioner, which of course is not recommended for general use.
667: Input Parameters:
668: + pc - preconditioner
669: - x - input vector
671: Output Parameter:
672: . y - preconditioned vector
673: */
674: PetscErrorCode MatrixFreePreconditioner(PC pc,Vec x,Vec y)
675: {
677: VecCopy(x,y);
678: return 0;
679: }
682: /*TEST
684: test:
685: args: -nox -snes_monitor_cancel -snes_monitor_short -ksp_gmres_cgs_refinement_type refine_always
687: test:
688: suffix: 2
689: nsize: 3
690: args: -nox -pc_type asm -mat_type mpiaij -snes_monitor_cancel -snes_monitor_short -ksp_gmres_cgs_refinement_type refine_always
692: test:
693: suffix: 3
694: nsize: 2
695: args: -nox -snes_monitor_cancel -snes_monitor_short -ksp_gmres_cgs_refinement_type refine_always
697: test:
698: suffix: 4
699: args: -nox -pre_check_iterates -post_check_iterates
701: TEST*/