Actual source code: ex3.c
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: /*
11: Include "petscdm.h" so that we can use data management objects (DMs)
12: Include "petscdmda.h" so that we can use distributed arrays (DMDAs).
13: Include "petscsnes.h" so that we can use SNES solvers. Note that this
14: file automatically includes:
15: petscsys.h - base PETSc routines
16: petscvec.h - vectors
17: petscmat.h - matrices
18: petscis.h - index sets
19: petscksp.h - Krylov subspace methods
20: petscviewer.h - viewers
21: petscpc.h - preconditioners
22: petscksp.h - linear solvers
23: */
25: #include <petscdm.h>
26: #include <petscdmda.h>
27: #include <petscsnes.h>
29: /*
30: User-defined routines.
31: */
32: PetscErrorCode FormJacobian(SNES,Vec,Mat,Mat,void*);
33: PetscErrorCode FormFunction(SNES,Vec,Vec,void*);
34: PetscErrorCode FormInitialGuess(Vec);
35: PetscErrorCode Monitor(SNES,PetscInt,PetscReal,void*);
36: PetscErrorCode PreCheck(SNESLineSearch,Vec,Vec,PetscBool*,void*);
37: PetscErrorCode PostCheck(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*);
38: PetscErrorCode PostSetSubKSP(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*);
39: PetscErrorCode MatrixFreePreconditioner(PC,Vec,Vec);
41: /*
42: User-defined application context
43: */
44: typedef struct {
45: DM da; /* distributed array */
46: Vec F; /* right-hand-side of PDE */
47: PetscMPIInt rank; /* rank of processor */
48: PetscMPIInt size; /* size of communicator */
49: PetscReal h; /* mesh spacing */
50: PetscBool sjerr; /* if or not to test jacobian domain error */
51: } ApplicationCtx;
53: /*
54: User-defined context for monitoring
55: */
56: typedef struct {
57: PetscViewer viewer;
58: } MonitorCtx;
60: /*
61: User-defined context for checking candidate iterates that are
62: determined by line search methods
63: */
64: typedef struct {
65: Vec last_step; /* previous iterate */
66: PetscReal tolerance; /* tolerance for changes between successive iterates */
67: ApplicationCtx *user;
68: } StepCheckCtx;
70: typedef struct {
71: PetscInt its0; /* num of prevous outer KSP iterations */
72: } SetSubKSPCtx;
74: int main(int argc,char **argv)
75: {
76: SNES snes; /* SNES context */
77: SNESLineSearch linesearch; /* SNESLineSearch context */
78: Mat J; /* Jacobian matrix */
79: ApplicationCtx ctx; /* user-defined context */
80: Vec x,r,U,F; /* vectors */
81: MonitorCtx monP; /* monitoring context */
82: StepCheckCtx checkP; /* step-checking context */
83: SetSubKSPCtx checkP1;
84: PetscBool pre_check,post_check,post_setsubksp; /* flag indicating whether we're checking candidate iterates */
85: PetscScalar xp,*FF,*UU,none = -1.0;
86: PetscInt its,N = 5,i,maxit,maxf,xs,xm;
87: PetscReal abstol,rtol,stol,norm;
88: PetscBool flg,viewinitial = PETSC_FALSE;
90: PetscInitialize(&argc,&argv,(char*)0,help);
91: MPI_Comm_rank(PETSC_COMM_WORLD,&ctx.rank);
92: MPI_Comm_size(PETSC_COMM_WORLD,&ctx.size);
93: PetscOptionsGetInt(NULL,NULL,"-n",&N,NULL);
94: ctx.h = 1.0/(N-1);
95: ctx.sjerr = PETSC_FALSE;
96: PetscOptionsGetBool(NULL,NULL,"-test_jacobian_domain_error",&ctx.sjerr,NULL);
97: PetscOptionsGetBool(NULL,NULL,"-view_initial",&viewinitial,NULL);
99: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
100: Create nonlinear solver context
101: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
103: SNESCreate(PETSC_COMM_WORLD,&snes);
105: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
106: Create vector data structures; set function evaluation routine
107: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
109: /*
110: Create distributed array (DMDA) to manage parallel grid and vectors
111: */
112: DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,N,1,1,NULL,&ctx.da);
113: DMSetFromOptions(ctx.da);
114: DMSetUp(ctx.da);
116: /*
117: Extract global and local vectors from DMDA; then duplicate for remaining
118: vectors that are the same types
119: */
120: DMCreateGlobalVector(ctx.da,&x);
121: VecDuplicate(x,&r);
122: VecDuplicate(x,&F); ctx.F = F;
123: VecDuplicate(x,&U);
125: /*
126: Set function evaluation routine and vector. Whenever the nonlinear
127: solver needs to compute the nonlinear function, it will call this
128: routine.
129: - Note that the final routine argument is the user-defined
130: context that provides application-specific data for the
131: function evaluation routine.
132: */
133: SNESSetFunction(snes,r,FormFunction,&ctx);
135: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
136: Create matrix data structure; set Jacobian evaluation routine
137: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
139: MatCreate(PETSC_COMM_WORLD,&J);
140: MatSetSizes(J,PETSC_DECIDE,PETSC_DECIDE,N,N);
141: MatSetFromOptions(J);
142: MatSeqAIJSetPreallocation(J,3,NULL);
143: MatMPIAIJSetPreallocation(J,3,NULL,3,NULL);
145: /*
146: Set Jacobian matrix data structure and default Jacobian evaluation
147: routine. Whenever the nonlinear solver needs to compute the
148: Jacobian matrix, it will call this routine.
149: - Note that the final routine argument is the user-defined
150: context that provides application-specific data for the
151: Jacobian evaluation routine.
152: */
153: SNESSetJacobian(snes,J,J,FormJacobian,&ctx);
155: /*
156: Optionally allow user-provided preconditioner
157: */
158: PetscOptionsHasName(NULL,NULL,"-user_precond",&flg);
159: if (flg) {
160: KSP ksp;
161: PC pc;
162: SNESGetKSP(snes,&ksp);
163: KSPGetPC(ksp,&pc);
164: PCSetType(pc,PCSHELL);
165: PCShellSetApply(pc,MatrixFreePreconditioner);
166: }
168: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
169: Customize nonlinear solver; set runtime options
170: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
172: /*
173: Set an optional user-defined monitoring routine
174: */
175: PetscViewerDrawOpen(PETSC_COMM_WORLD,0,0,0,0,400,400,&monP.viewer);
176: SNESMonitorSet(snes,Monitor,&monP,0);
178: /*
179: Set names for some vectors to facilitate monitoring (optional)
180: */
181: PetscObjectSetName((PetscObject)x,"Approximate Solution");
182: PetscObjectSetName((PetscObject)U,"Exact Solution");
184: /*
185: Set SNES/KSP/KSP/PC runtime options, e.g.,
186: -snes_view -snes_monitor -ksp_type <ksp> -pc_type <pc>
187: */
188: SNESSetFromOptions(snes);
190: /*
191: Set an optional user-defined routine to check the validity of candidate
192: iterates that are determined by line search methods
193: */
194: SNESGetLineSearch(snes, &linesearch);
195: PetscOptionsHasName(NULL,NULL,"-post_check_iterates",&post_check);
197: if (post_check) {
198: PetscPrintf(PETSC_COMM_WORLD,"Activating post step checking routine\n");
199: SNESLineSearchSetPostCheck(linesearch,PostCheck,&checkP);
200: VecDuplicate(x,&(checkP.last_step));
202: checkP.tolerance = 1.0;
203: checkP.user = &ctx;
205: PetscOptionsGetReal(NULL,NULL,"-check_tol",&checkP.tolerance,NULL);
206: }
208: PetscOptionsHasName(NULL,NULL,"-post_setsubksp",&post_setsubksp);
209: if (post_setsubksp) {
210: PetscPrintf(PETSC_COMM_WORLD,"Activating post setsubksp\n");
211: SNESLineSearchSetPostCheck(linesearch,PostSetSubKSP,&checkP1);
212: }
214: PetscOptionsHasName(NULL,NULL,"-pre_check_iterates",&pre_check);
215: if (pre_check) {
216: PetscPrintf(PETSC_COMM_WORLD,"Activating pre step checking routine\n");
217: SNESLineSearchSetPreCheck(linesearch,PreCheck,&checkP);
218: }
220: /*
221: Print parameters used for convergence testing (optional) ... just
222: to demonstrate this routine; this information is also printed with
223: the option -snes_view
224: */
225: SNESGetTolerances(snes,&abstol,&rtol,&stol,&maxit,&maxf);
226: PetscPrintf(PETSC_COMM_WORLD,"atol=%g, rtol=%g, stol=%g, maxit=%D, maxf=%D\n",(double)abstol,(double)rtol,(double)stol,maxit,maxf);
228: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
229: Initialize application:
230: Store right-hand-side of PDE and exact solution
231: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
233: /*
234: Get local grid boundaries (for 1-dimensional DMDA):
235: xs, xm - starting grid index, width of local grid (no ghost points)
236: */
237: DMDAGetCorners(ctx.da,&xs,NULL,NULL,&xm,NULL,NULL);
239: /*
240: Get pointers to vector data
241: */
242: DMDAVecGetArray(ctx.da,F,&FF);
243: DMDAVecGetArray(ctx.da,U,&UU);
245: /*
246: Compute local vector entries
247: */
248: xp = ctx.h*xs;
249: for (i=xs; i<xs+xm; i++) {
250: FF[i] = 6.0*xp + PetscPowScalar(xp+1.e-12,6.0); /* +1.e-12 is to prevent 0^6 */
251: UU[i] = xp*xp*xp;
252: xp += ctx.h;
253: }
255: /*
256: Restore vectors
257: */
258: DMDAVecRestoreArray(ctx.da,F,&FF);
259: DMDAVecRestoreArray(ctx.da,U,&UU);
260: if (viewinitial) {
261: VecView(U,0);
262: VecView(F,0);
263: }
265: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
266: Evaluate initial guess; then solve nonlinear system
267: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
269: /*
270: Note: The user should initialize the vector, x, with the initial guess
271: for the nonlinear solver prior to calling SNESSolve(). In particular,
272: to employ an initial guess of zero, the user should explicitly set
273: this vector to zero by calling VecSet().
274: */
275: FormInitialGuess(x);
276: SNESSolve(snes,NULL,x);
277: SNESGetIterationNumber(snes,&its);
278: PetscPrintf(PETSC_COMM_WORLD,"Number of SNES iterations = %D\n",its);
280: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
281: Check solution and clean up
282: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
284: /*
285: Check the error
286: */
287: VecAXPY(x,none,U);
288: VecNorm(x,NORM_2,&norm);
289: PetscPrintf(PETSC_COMM_WORLD,"Norm of error %g Iterations %D\n",(double)norm,its);
290: if (ctx.sjerr) {
291: SNESType snestype;
292: SNESGetType(snes,&snestype);
293: PetscPrintf(PETSC_COMM_WORLD,"SNES Type %s\n",snestype);
294: }
296: /*
297: Free work space. All PETSc objects should be destroyed when they
298: are no longer needed.
299: */
300: PetscViewerDestroy(&monP.viewer);
301: if (post_check) VecDestroy(&checkP.last_step);
302: VecDestroy(&x);
303: VecDestroy(&r);
304: VecDestroy(&U);
305: VecDestroy(&F);
306: MatDestroy(&J);
307: SNESDestroy(&snes);
308: DMDestroy(&ctx.da);
309: PetscFinalize();
310: return 0;
311: }
313: /* ------------------------------------------------------------------- */
314: /*
315: FormInitialGuess - Computes initial guess.
317: Input/Output Parameter:
318: . x - the solution vector
319: */
320: PetscErrorCode FormInitialGuess(Vec x)
321: {
322: PetscScalar pfive = .50;
325: VecSet(x,pfive);
326: return 0;
327: }
329: /* ------------------------------------------------------------------- */
330: /*
331: FormFunction - Evaluates nonlinear function, F(x).
333: Input Parameters:
334: . snes - the SNES context
335: . x - input vector
336: . ctx - optional user-defined context, as set by SNESSetFunction()
338: Output Parameter:
339: . f - function vector
341: Note:
342: The user-defined context can contain any application-specific
343: data needed for the function evaluation.
344: */
345: PetscErrorCode FormFunction(SNES snes,Vec x,Vec f,void *ctx)
346: {
347: ApplicationCtx *user = (ApplicationCtx*) ctx;
348: DM da = user->da;
349: PetscScalar *xx,*ff,*FF,d;
350: PetscInt i,M,xs,xm;
351: Vec xlocal;
354: DMGetLocalVector(da,&xlocal);
355: /*
356: Scatter ghost points to local vector, using the 2-step process
357: DMGlobalToLocalBegin(), DMGlobalToLocalEnd().
358: By placing code between these two statements, computations can
359: be done while messages are in transition.
360: */
361: DMGlobalToLocalBegin(da,x,INSERT_VALUES,xlocal);
362: DMGlobalToLocalEnd(da,x,INSERT_VALUES,xlocal);
364: /*
365: Get pointers to vector data.
366: - The vector xlocal includes ghost point; the vectors x and f do
367: NOT include ghost points.
368: - Using DMDAVecGetArray() allows accessing the values using global ordering
369: */
370: DMDAVecGetArray(da,xlocal,&xx);
371: DMDAVecGetArray(da,f,&ff);
372: DMDAVecGetArray(da,user->F,&FF);
374: /*
375: Get local grid boundaries (for 1-dimensional DMDA):
376: xs, xm - starting grid index, width of local grid (no ghost points)
377: */
378: DMDAGetCorners(da,&xs,NULL,NULL,&xm,NULL,NULL);
379: DMDAGetInfo(da,NULL,&M,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);
381: /*
382: Set function values for boundary points; define local interior grid point range:
383: xsi - starting interior grid index
384: xei - ending interior grid index
385: */
386: if (xs == 0) { /* left boundary */
387: ff[0] = xx[0];
388: xs++;xm--;
389: }
390: if (xs+xm == M) { /* right boundary */
391: ff[xs+xm-1] = xx[xs+xm-1] - 1.0;
392: xm--;
393: }
395: /*
396: Compute function over locally owned part of the grid (interior points only)
397: */
398: d = 1.0/(user->h*user->h);
399: 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];
401: /*
402: Restore vectors
403: */
404: DMDAVecRestoreArray(da,xlocal,&xx);
405: DMDAVecRestoreArray(da,f,&ff);
406: DMDAVecRestoreArray(da,user->F,&FF);
407: DMRestoreLocalVector(da,&xlocal);
408: return 0;
409: }
411: /* ------------------------------------------------------------------- */
412: /*
413: FormJacobian - Evaluates Jacobian matrix.
415: Input Parameters:
416: . snes - the SNES context
417: . x - input vector
418: . dummy - optional user-defined context (not used here)
420: Output Parameters:
421: . jac - Jacobian matrix
422: . B - optionally different preconditioning matrix
423: . flag - flag indicating matrix structure
424: */
425: PetscErrorCode FormJacobian(SNES snes,Vec x,Mat jac,Mat B,void *ctx)
426: {
427: ApplicationCtx *user = (ApplicationCtx*) ctx;
428: PetscScalar *xx,d,A[3];
429: PetscInt i,j[3],M,xs,xm;
430: DM da = user->da;
433: if (user->sjerr) {
434: SNESSetJacobianDomainError(snes);
435: return 0;
436: }
437: /*
438: Get pointer to vector data
439: */
440: DMDAVecGetArrayRead(da,x,&xx);
441: DMDAGetCorners(da,&xs,NULL,NULL,&xm,NULL,NULL);
443: /*
444: Get range of locally owned matrix
445: */
446: DMDAGetInfo(da,NULL,&M,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);
448: /*
449: Determine starting and ending local indices for interior grid points.
450: Set Jacobian entries for boundary points.
451: */
453: if (xs == 0) { /* left boundary */
454: i = 0; A[0] = 1.0;
456: MatSetValues(jac,1,&i,1,&i,A,INSERT_VALUES);
457: xs++;xm--;
458: }
459: if (xs+xm == M) { /* right boundary */
460: i = M-1;
461: A[0] = 1.0;
462: MatSetValues(jac,1,&i,1,&i,A,INSERT_VALUES);
463: xm--;
464: }
466: /*
467: Interior grid points
468: - Note that in this case we set all elements for a particular
469: row at once.
470: */
471: d = 1.0/(user->h*user->h);
472: for (i=xs; i<xs+xm; i++) {
473: j[0] = i - 1; j[1] = i; j[2] = i + 1;
474: A[0] = A[2] = d; A[1] = -2.0*d + 2.0*xx[i];
475: MatSetValues(jac,1,&i,3,j,A,INSERT_VALUES);
476: }
478: /*
479: Assemble matrix, using the 2-step process:
480: MatAssemblyBegin(), MatAssemblyEnd().
481: By placing code between these two statements, computations can be
482: done while messages are in transition.
484: Also, restore vector.
485: */
487: MatAssemblyBegin(jac,MAT_FINAL_ASSEMBLY);
488: DMDAVecRestoreArrayRead(da,x,&xx);
489: MatAssemblyEnd(jac,MAT_FINAL_ASSEMBLY);
491: return 0;
492: }
494: /* ------------------------------------------------------------------- */
495: /*
496: Monitor - Optional user-defined monitoring routine that views the
497: current iterate with an x-window plot. Set by SNESMonitorSet().
499: Input Parameters:
500: snes - the SNES context
501: its - iteration number
502: norm - 2-norm function value (may be estimated)
503: ctx - optional user-defined context for private data for the
504: monitor routine, as set by SNESMonitorSet()
506: Note:
507: See the manpage for PetscViewerDrawOpen() for useful runtime options,
508: such as -nox to deactivate all x-window output.
509: */
510: PetscErrorCode Monitor(SNES snes,PetscInt its,PetscReal fnorm,void *ctx)
511: {
512: MonitorCtx *monP = (MonitorCtx*) ctx;
513: Vec x;
516: PetscPrintf(PETSC_COMM_WORLD,"iter = %D,SNES Function norm %g\n",its,(double)fnorm);
517: SNESGetSolution(snes,&x);
518: VecView(x,monP->viewer);
519: return 0;
520: }
522: /* ------------------------------------------------------------------- */
523: /*
524: PreCheck - Optional user-defined routine that checks the validity of
525: candidate steps of a line search method. Set by SNESLineSearchSetPreCheck().
527: Input Parameters:
528: snes - the SNES context
529: xcurrent - current solution
530: y - search direction and length
532: Output Parameters:
533: y - proposed step (search direction and length) (possibly changed)
534: changed_y - tells if the step has changed or not
535: */
536: PetscErrorCode PreCheck(SNESLineSearch linesearch,Vec xcurrent,Vec y, PetscBool *changed_y, void * ctx)
537: {
539: *changed_y = PETSC_FALSE;
540: return 0;
541: }
543: /* ------------------------------------------------------------------- */
544: /*
545: PostCheck - Optional user-defined routine that checks the validity of
546: candidate steps of a line search method. Set by SNESLineSearchSetPostCheck().
548: Input Parameters:
549: snes - the SNES context
550: ctx - optional user-defined context for private data for the
551: monitor routine, as set by SNESLineSearchSetPostCheck()
552: xcurrent - current solution
553: y - search direction and length
554: x - the new candidate iterate
556: Output Parameters:
557: y - proposed step (search direction and length) (possibly changed)
558: x - current iterate (possibly modified)
560: */
561: PetscErrorCode PostCheck(SNESLineSearch linesearch,Vec xcurrent,Vec y,Vec x,PetscBool *changed_y,PetscBool *changed_x, void * ctx)
562: {
563: PetscInt i,iter,xs,xm;
564: StepCheckCtx *check;
565: ApplicationCtx *user;
566: PetscScalar *xa,*xa_last,tmp;
567: PetscReal rdiff;
568: DM da;
569: SNES snes;
572: *changed_x = PETSC_FALSE;
573: *changed_y = PETSC_FALSE;
575: SNESLineSearchGetSNES(linesearch, &snes);
576: check = (StepCheckCtx*)ctx;
577: user = check->user;
578: SNESGetIterationNumber(snes,&iter);
580: /* iteration 1 indicates we are working on the second iteration */
581: if (iter > 0) {
582: da = user->da;
583: PetscPrintf(PETSC_COMM_WORLD,"Checking candidate step at iteration %D with tolerance %g\n",iter,(double)check->tolerance);
585: /* Access local array data */
586: DMDAVecGetArray(da,check->last_step,&xa_last);
587: DMDAVecGetArray(da,x,&xa);
588: DMDAGetCorners(da,&xs,NULL,NULL,&xm,NULL,NULL);
590: /*
591: If we fail the user-defined check for validity of the candidate iterate,
592: then modify the iterate as we like. (Note that the particular modification
593: below is intended simply to demonstrate how to manipulate this data, not
594: as a meaningful or appropriate choice.)
595: */
596: for (i=xs; i<xs+xm; i++) {
597: if (!PetscAbsScalar(xa[i])) rdiff = 2*check->tolerance;
598: else rdiff = PetscAbsScalar((xa[i] - xa_last[i])/xa[i]);
599: if (rdiff > check->tolerance) {
600: tmp = xa[i];
601: xa[i] = .5*(xa[i] + xa_last[i]);
602: *changed_x = PETSC_TRUE;
603: 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]));
604: }
605: }
606: DMDAVecRestoreArray(da,check->last_step,&xa_last);
607: DMDAVecRestoreArray(da,x,&xa);
608: }
609: VecCopy(x,check->last_step);
610: return 0;
611: }
613: /* ------------------------------------------------------------------- */
614: /*
615: PostSetSubKSP - Optional user-defined routine that reset SubKSP options when hierarchical bjacobi PC is used
616: e.g,
617: 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
618: Set by SNESLineSearchSetPostCheck().
620: Input Parameters:
621: linesearch - the LineSearch context
622: xcurrent - current solution
623: y - search direction and length
624: x - the new candidate iterate
626: Output Parameters:
627: y - proposed step (search direction and length) (possibly changed)
628: x - current iterate (possibly modified)
630: */
631: PetscErrorCode PostSetSubKSP(SNESLineSearch linesearch,Vec xcurrent,Vec y,Vec x,PetscBool *changed_y,PetscBool *changed_x, void * ctx)
632: {
633: SetSubKSPCtx *check;
634: PetscInt iter,its,sub_its,maxit;
635: KSP ksp,sub_ksp,*sub_ksps;
636: PC pc;
637: PetscReal ksp_ratio;
638: SNES snes;
641: SNESLineSearchGetSNES(linesearch, &snes);
642: check = (SetSubKSPCtx*)ctx;
643: SNESGetIterationNumber(snes,&iter);
644: SNESGetKSP(snes,&ksp);
645: KSPGetPC(ksp,&pc);
646: PCBJacobiGetSubKSP(pc,NULL,NULL,&sub_ksps);
647: sub_ksp = sub_ksps[0];
648: KSPGetIterationNumber(ksp,&its); /* outer KSP iteration number */
649: KSPGetIterationNumber(sub_ksp,&sub_its); /* inner KSP iteration number */
651: if (iter) {
652: PetscPrintf(PETSC_COMM_WORLD," ...PostCheck snes iteration %D, ksp_it %D %D, subksp_it %D\n",iter,check->its0,its,sub_its);
653: ksp_ratio = ((PetscReal)(its))/check->its0;
654: maxit = (PetscInt)(ksp_ratio*sub_its + 0.5);
655: if (maxit < 2) maxit = 2;
656: KSPSetTolerances(sub_ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,maxit);
657: PetscPrintf(PETSC_COMM_WORLD," ...ksp_ratio %g, new maxit %D\n\n",(double)ksp_ratio,maxit);
658: }
659: check->its0 = its; /* save current outer KSP iteration number */
660: return 0;
661: }
663: /* ------------------------------------------------------------------- */
664: /*
665: MatrixFreePreconditioner - This routine demonstrates the use of a
666: user-provided preconditioner. This code implements just the null
667: preconditioner, which of course is not recommended for general use.
669: Input Parameters:
670: + pc - preconditioner
671: - x - input vector
673: Output Parameter:
674: . y - preconditioned vector
675: */
676: PetscErrorCode MatrixFreePreconditioner(PC pc,Vec x,Vec y)
677: {
678: VecCopy(x,y);
679: return 0;
680: }
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:
702: suffix: 5
703: requires: double !complex !single
704: nsize: 2
705: args: -nox -snes_test_jacobian -snes_test_jacobian_view
707: test:
708: suffix: 6
709: requires: double !complex !single
710: nsize: 4
711: args: -test_jacobian_domain_error -snes_converged_reason -snes_check_jacobian_domain_error 1
713: test:
714: suffix: 7
715: requires: double !complex !single
716: nsize: 4
717: args: -test_jacobian_domain_error -snes_converged_reason -snes_type newtontr -snes_check_jacobian_domain_error 1
719: test:
720: suffix: 8
721: requires: double !complex !single
722: nsize: 4
723: args: -test_jacobian_domain_error -snes_converged_reason -snes_type vinewtonrsls -snes_check_jacobian_domain_error 1
725: test:
726: suffix: 9
727: requires: double !complex !single
728: nsize: 4
729: args: -test_jacobian_domain_error -snes_converged_reason -snes_type vinewtonssls -snes_check_jacobian_domain_error 1
731: test:
732: suffix: 10
733: requires: double !complex !single
734: nsize: 4
735: args: -test_jacobian_domain_error -snes_converged_reason -snes_type qn -snes_qn_scale_type jacobian -snes_check_jacobian_domain_error 1
737: test:
738: suffix: 11
739: requires: double !complex !single
740: nsize: 4
741: args: -test_jacobian_domain_error -snes_converged_reason -snes_type ms -snes_ms_type m62 -snes_ms_damping 0.9 -snes_check_jacobian_domain_error 1
743: test:
744: suffix: 12
745: args: -view_initial
746: filter: grep -v "type:"
748: test:
749: suffix: 13
750: requires: double !complex !single
751: nsize: 4
752: args: -test_jacobian_domain_error -snes_converged_reason -snes_type newtontrdc -snes_check_jacobian_domain_error 1
754: TEST*/