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