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;
93: PetscInt its,N = 5,i,maxit,maxf,xs,xm;
94: PetscReal abstol,rtol,stol,norm;
95: PetscBool flg,viewinitial = PETSC_FALSE;
97: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
98: MPI_Comm_rank(PETSC_COMM_WORLD,&ctx.rank);
99: MPI_Comm_size(PETSC_COMM_WORLD,&ctx.size);
100: PetscOptionsGetInt(NULL,NULL,"-n",&N,NULL);
101: ctx.h = 1.0/(N-1);
102: ctx.sjerr = PETSC_FALSE;
103: PetscOptionsGetBool(NULL,NULL,"-test_jacobian_domain_error",&ctx.sjerr,NULL);
104: PetscOptionsGetBool(NULL,NULL,"-view_initial",&viewinitial,NULL);
106: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
107: Create nonlinear solver context
108: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
110: SNESCreate(PETSC_COMM_WORLD,&snes);
112: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
113: Create vector data structures; set function evaluation routine
114: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
116: /*
117: Create distributed array (DMDA) to manage parallel grid and vectors
118: */
119: DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,N,1,1,NULL,&ctx.da);
120: DMSetFromOptions(ctx.da);
121: DMSetUp(ctx.da);
123: /*
124: Extract global and local vectors from DMDA; then duplicate for remaining
125: vectors that are the same types
126: */
127: DMCreateGlobalVector(ctx.da,&x);
128: VecDuplicate(x,&r);
129: VecDuplicate(x,&F); ctx.F = F;
130: VecDuplicate(x,&U);
132: /*
133: Set function evaluation routine and vector. Whenever the nonlinear
134: solver needs to compute the nonlinear function, it will call this
135: routine.
136: - Note that the final routine argument is the user-defined
137: context that provides application-specific data for the
138: function evaluation routine.
139: */
140: SNESSetFunction(snes,r,FormFunction,&ctx);
142: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
143: Create matrix data structure; set Jacobian evaluation routine
144: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
146: MatCreate(PETSC_COMM_WORLD,&J);
147: MatSetSizes(J,PETSC_DECIDE,PETSC_DECIDE,N,N);
148: MatSetFromOptions(J);
149: MatSeqAIJSetPreallocation(J,3,NULL);
150: MatMPIAIJSetPreallocation(J,3,NULL,3,NULL);
152: /*
153: Set Jacobian matrix data structure and default Jacobian evaluation
154: routine. Whenever the nonlinear solver needs to compute the
155: Jacobian matrix, it will call this routine.
156: - Note that the final routine argument is the user-defined
157: context that provides application-specific data for the
158: Jacobian evaluation routine.
159: */
160: SNESSetJacobian(snes,J,J,FormJacobian,&ctx);
162: /*
163: Optionally allow user-provided preconditioner
164: */
165: PetscOptionsHasName(NULL,NULL,"-user_precond",&flg);
166: if (flg) {
167: KSP ksp;
168: PC pc;
169: SNESGetKSP(snes,&ksp);
170: KSPGetPC(ksp,&pc);
171: PCSetType(pc,PCSHELL);
172: PCShellSetApply(pc,MatrixFreePreconditioner);
173: }
175: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
176: Customize nonlinear solver; set runtime options
177: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
179: /*
180: Set an optional user-defined monitoring routine
181: */
182: PetscViewerDrawOpen(PETSC_COMM_WORLD,0,0,0,0,400,400,&monP.viewer);
183: SNESMonitorSet(snes,Monitor,&monP,0);
185: /*
186: Set names for some vectors to facilitate monitoring (optional)
187: */
188: PetscObjectSetName((PetscObject)x,"Approximate Solution");
189: PetscObjectSetName((PetscObject)U,"Exact Solution");
191: /*
192: Set SNES/KSP/KSP/PC runtime options, e.g.,
193: -snes_view -snes_monitor -ksp_type <ksp> -pc_type <pc>
194: */
195: SNESSetFromOptions(snes);
197: /*
198: Set an optional user-defined routine to check the validity of candidate
199: iterates that are determined by line search methods
200: */
201: SNESGetLineSearch(snes, &linesearch);
202: PetscOptionsHasName(NULL,NULL,"-post_check_iterates",&post_check);
204: if (post_check) {
205: PetscPrintf(PETSC_COMM_WORLD,"Activating post step checking routine\n");
206: SNESLineSearchSetPostCheck(linesearch,PostCheck,&checkP);
207: VecDuplicate(x,&(checkP.last_step));
209: checkP.tolerance = 1.0;
210: checkP.user = &ctx;
212: PetscOptionsGetReal(NULL,NULL,"-check_tol",&checkP.tolerance,NULL);
213: }
215: PetscOptionsHasName(NULL,NULL,"-post_setsubksp",&post_setsubksp);
216: if (post_setsubksp) {
217: PetscPrintf(PETSC_COMM_WORLD,"Activating post setsubksp\n");
218: SNESLineSearchSetPostCheck(linesearch,PostSetSubKSP,&checkP1);
219: }
221: PetscOptionsHasName(NULL,NULL,"-pre_check_iterates",&pre_check);
222: if (pre_check) {
223: PetscPrintf(PETSC_COMM_WORLD,"Activating pre step checking routine\n");
224: SNESLineSearchSetPreCheck(linesearch,PreCheck,&checkP);
225: }
227: /*
228: Print parameters used for convergence testing (optional) ... just
229: to demonstrate this routine; this information is also printed with
230: the option -snes_view
231: */
232: SNESGetTolerances(snes,&abstol,&rtol,&stol,&maxit,&maxf);
233: PetscPrintf(PETSC_COMM_WORLD,"atol=%g, rtol=%g, stol=%g, maxit=%D, maxf=%D\n",(double)abstol,(double)rtol,(double)stol,maxit,maxf);
235: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
236: Initialize application:
237: Store right-hand-side of PDE and exact solution
238: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
240: /*
241: Get local grid boundaries (for 1-dimensional DMDA):
242: xs, xm - starting grid index, width of local grid (no ghost points)
243: */
244: DMDAGetCorners(ctx.da,&xs,NULL,NULL,&xm,NULL,NULL);
246: /*
247: Get pointers to vector data
248: */
249: DMDAVecGetArray(ctx.da,F,&FF);
250: DMDAVecGetArray(ctx.da,U,&UU);
252: /*
253: Compute local vector entries
254: */
255: xp = ctx.h*xs;
256: for (i=xs; i<xs+xm; i++) {
257: FF[i] = 6.0*xp + PetscPowScalar(xp+1.e-12,6.0); /* +1.e-12 is to prevent 0^6 */
258: UU[i] = xp*xp*xp;
259: xp += ctx.h;
260: }
262: /*
263: Restore vectors
264: */
265: DMDAVecRestoreArray(ctx.da,F,&FF);
266: DMDAVecRestoreArray(ctx.da,U,&UU);
267: if (viewinitial) {
268: VecView(U,0);
269: VecView(F,0);
270: }
272: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
273: Evaluate initial guess; then solve nonlinear system
274: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
276: /*
277: Note: The user should initialize the vector, x, with the initial guess
278: for the nonlinear solver prior to calling SNESSolve(). In particular,
279: to employ an initial guess of zero, the user should explicitly set
280: this vector to zero by calling VecSet().
281: */
282: FormInitialGuess(x);
283: SNESSolve(snes,NULL,x);
284: SNESGetIterationNumber(snes,&its);
285: PetscPrintf(PETSC_COMM_WORLD,"Number of SNES iterations = %D\n",its);
287: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
288: Check solution and clean up
289: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
291: /*
292: Check the error
293: */
294: VecAXPY(x,none,U);
295: VecNorm(x,NORM_2,&norm);
296: PetscPrintf(PETSC_COMM_WORLD,"Norm of error %g Iterations %D\n",(double)norm,its);
297: if (ctx.sjerr) {
298: SNESType snestype;
299: SNESGetType(snes,&snestype);
300: PetscPrintf(PETSC_COMM_WORLD,"SNES Type %s\n",snestype);
301: }
303: /*
304: Free work space. All PETSc objects should be destroyed when they
305: are no longer needed.
306: */
307: PetscViewerDestroy(&monP.viewer);
308: if (post_check) {VecDestroy(&checkP.last_step);}
309: VecDestroy(&x);
310: VecDestroy(&r);
311: VecDestroy(&U);
312: VecDestroy(&F);
313: MatDestroy(&J);
314: SNESDestroy(&snes);
315: DMDestroy(&ctx.da);
316: PetscFinalize();
317: return ierr;
318: }
320: /* ------------------------------------------------------------------- */
321: /*
322: FormInitialGuess - Computes initial guess.
324: Input/Output Parameter:
325: . x - the solution vector
326: */
327: PetscErrorCode FormInitialGuess(Vec x)
328: {
330: PetscScalar pfive = .50;
333: VecSet(x,pfive);
334: return(0);
335: }
337: /* ------------------------------------------------------------------- */
338: /*
339: FormFunction - Evaluates nonlinear function, F(x).
341: Input Parameters:
342: . snes - the SNES context
343: . x - input vector
344: . ctx - optional user-defined context, as set by SNESSetFunction()
346: Output Parameter:
347: . f - function vector
349: Note:
350: The user-defined context can contain any application-specific
351: data needed for the function evaluation.
352: */
353: PetscErrorCode FormFunction(SNES snes,Vec x,Vec f,void *ctx)
354: {
355: ApplicationCtx *user = (ApplicationCtx*) ctx;
356: DM da = user->da;
357: PetscScalar *xx,*ff,*FF,d;
359: PetscInt i,M,xs,xm;
360: Vec xlocal;
363: DMGetLocalVector(da,&xlocal);
364: /*
365: Scatter ghost points to local vector, using the 2-step process
366: DMGlobalToLocalBegin(), DMGlobalToLocalEnd().
367: By placing code between these two statements, computations can
368: be done while messages are in transition.
369: */
370: DMGlobalToLocalBegin(da,x,INSERT_VALUES,xlocal);
371: DMGlobalToLocalEnd(da,x,INSERT_VALUES,xlocal);
373: /*
374: Get pointers to vector data.
375: - The vector xlocal includes ghost point; the vectors x and f do
376: NOT include ghost points.
377: - Using DMDAVecGetArray() allows accessing the values using global ordering
378: */
379: DMDAVecGetArray(da,xlocal,&xx);
380: DMDAVecGetArray(da,f,&ff);
381: DMDAVecGetArray(da,user->F,&FF);
383: /*
384: Get local grid boundaries (for 1-dimensional DMDA):
385: xs, xm - starting grid index, width of local grid (no ghost points)
386: */
387: DMDAGetCorners(da,&xs,NULL,NULL,&xm,NULL,NULL);
388: DMDAGetInfo(da,NULL,&M,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);
390: /*
391: Set function values for boundary points; define local interior grid point range:
392: xsi - starting interior grid index
393: xei - ending interior grid index
394: */
395: if (xs == 0) { /* left boundary */
396: ff[0] = xx[0];
397: xs++;xm--;
398: }
399: if (xs+xm == M) { /* right boundary */
400: ff[xs+xm-1] = xx[xs+xm-1] - 1.0;
401: xm--;
402: }
404: /*
405: Compute function over locally owned part of the grid (interior points only)
406: */
407: d = 1.0/(user->h*user->h);
408: 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];
410: /*
411: Restore vectors
412: */
413: DMDAVecRestoreArray(da,xlocal,&xx);
414: DMDAVecRestoreArray(da,f,&ff);
415: DMDAVecRestoreArray(da,user->F,&FF);
416: DMRestoreLocalVector(da,&xlocal);
417: return(0);
418: }
420: /* ------------------------------------------------------------------- */
421: /*
422: FormJacobian - Evaluates Jacobian matrix.
424: Input Parameters:
425: . snes - the SNES context
426: . x - input vector
427: . dummy - optional user-defined context (not used here)
429: Output Parameters:
430: . jac - Jacobian matrix
431: . B - optionally different preconditioning matrix
432: . flag - flag indicating matrix structure
433: */
434: PetscErrorCode FormJacobian(SNES snes,Vec x,Mat jac,Mat B,void *ctx)
435: {
436: ApplicationCtx *user = (ApplicationCtx*) ctx;
437: PetscScalar *xx,d,A[3];
439: PetscInt i,j[3],M,xs,xm;
440: DM da = user->da;
443: if (user->sjerr) {
444: SNESSetJacobianDomainError(snes);
445: return(0);
446: }
447: /*
448: Get pointer to vector data
449: */
450: DMDAVecGetArrayRead(da,x,&xx);
451: DMDAGetCorners(da,&xs,NULL,NULL,&xm,NULL,NULL);
453: /*
454: Get range of locally owned matrix
455: */
456: DMDAGetInfo(da,NULL,&M,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL,NULL);
458: /*
459: Determine starting and ending local indices for interior grid points.
460: Set Jacobian entries for boundary points.
461: */
463: if (xs == 0) { /* left boundary */
464: i = 0; A[0] = 1.0;
466: MatSetValues(jac,1,&i,1,&i,A,INSERT_VALUES);
467: xs++;xm--;
468: }
469: if (xs+xm == M) { /* right boundary */
470: i = M-1;
471: A[0] = 1.0;
472: MatSetValues(jac,1,&i,1,&i,A,INSERT_VALUES);
473: xm--;
474: }
476: /*
477: Interior grid points
478: - Note that in this case we set all elements for a particular
479: row at once.
480: */
481: d = 1.0/(user->h*user->h);
482: for (i=xs; i<xs+xm; i++) {
483: j[0] = i - 1; j[1] = i; j[2] = i + 1;
484: A[0] = A[2] = d; A[1] = -2.0*d + 2.0*xx[i];
485: MatSetValues(jac,1,&i,3,j,A,INSERT_VALUES);
486: }
488: /*
489: Assemble matrix, using the 2-step process:
490: MatAssemblyBegin(), MatAssemblyEnd().
491: By placing code between these two statements, computations can be
492: done while messages are in transition.
494: Also, restore vector.
495: */
497: MatAssemblyBegin(jac,MAT_FINAL_ASSEMBLY);
498: DMDAVecRestoreArrayRead(da,x,&xx);
499: MatAssemblyEnd(jac,MAT_FINAL_ASSEMBLY);
501: return(0);
502: }
504: /* ------------------------------------------------------------------- */
505: /*
506: Monitor - Optional user-defined monitoring routine that views the
507: current iterate with an x-window plot. Set by SNESMonitorSet().
509: Input Parameters:
510: snes - the SNES context
511: its - iteration number
512: norm - 2-norm function value (may be estimated)
513: ctx - optional user-defined context for private data for the
514: monitor routine, as set by SNESMonitorSet()
516: Note:
517: See the manpage for PetscViewerDrawOpen() for useful runtime options,
518: such as -nox to deactivate all x-window output.
519: */
520: PetscErrorCode Monitor(SNES snes,PetscInt its,PetscReal fnorm,void *ctx)
521: {
523: MonitorCtx *monP = (MonitorCtx*) ctx;
524: Vec x;
527: PetscPrintf(PETSC_COMM_WORLD,"iter = %D,SNES Function norm %g\n",its,(double)fnorm);
528: SNESGetSolution(snes,&x);
529: VecView(x,monP->viewer);
530: return(0);
531: }
533: /* ------------------------------------------------------------------- */
534: /*
535: PreCheck - Optional user-defined routine that checks the validity of
536: candidate steps of a line search method. Set by SNESLineSearchSetPreCheck().
538: Input Parameters:
539: snes - the SNES context
540: xcurrent - current solution
541: y - search direction and length
543: Output Parameters:
544: y - proposed step (search direction and length) (possibly changed)
545: changed_y - tells if the step has changed or not
546: */
547: PetscErrorCode PreCheck(SNESLineSearch linesearch,Vec xcurrent,Vec y, PetscBool *changed_y, void * ctx)
548: {
550: *changed_y = PETSC_FALSE;
551: return(0);
552: }
554: /* ------------------------------------------------------------------- */
555: /*
556: PostCheck - Optional user-defined routine that checks the validity of
557: candidate steps of a line search method. Set by SNESLineSearchSetPostCheck().
559: Input Parameters:
560: snes - the SNES context
561: ctx - optional user-defined context for private data for the
562: monitor routine, as set by SNESLineSearchSetPostCheck()
563: xcurrent - current solution
564: y - search direction and length
565: x - the new candidate iterate
567: Output Parameters:
568: y - proposed step (search direction and length) (possibly changed)
569: x - current iterate (possibly modified)
571: */
572: PetscErrorCode PostCheck(SNESLineSearch linesearch,Vec xcurrent,Vec y,Vec x,PetscBool *changed_y,PetscBool *changed_x, void * ctx)
573: {
575: PetscInt i,iter,xs,xm;
576: StepCheckCtx *check;
577: ApplicationCtx *user;
578: PetscScalar *xa,*xa_last,tmp;
579: PetscReal rdiff;
580: DM da;
581: SNES snes;
584: *changed_x = PETSC_FALSE;
585: *changed_y = PETSC_FALSE;
587: SNESLineSearchGetSNES(linesearch, &snes);
588: check = (StepCheckCtx*)ctx;
589: user = check->user;
590: SNESGetIterationNumber(snes,&iter);
592: /* iteration 1 indicates we are working on the second iteration */
593: if (iter > 0) {
594: da = user->da;
595: PetscPrintf(PETSC_COMM_WORLD,"Checking candidate step at iteration %D with tolerance %g\n",iter,(double)check->tolerance);
597: /* Access local array data */
598: DMDAVecGetArray(da,check->last_step,&xa_last);
599: DMDAVecGetArray(da,x,&xa);
600: DMDAGetCorners(da,&xs,NULL,NULL,&xm,NULL,NULL);
602: /*
603: If we fail the user-defined check for validity of the candidate iterate,
604: then modify the iterate as we like. (Note that the particular modification
605: below is intended simply to demonstrate how to manipulate this data, not
606: as a meaningful or appropriate choice.)
607: */
608: for (i=xs; i<xs+xm; i++) {
609: if (!PetscAbsScalar(xa[i])) rdiff = 2*check->tolerance;
610: else rdiff = PetscAbsScalar((xa[i] - xa_last[i])/xa[i]);
611: if (rdiff > check->tolerance) {
612: tmp = xa[i];
613: xa[i] = .5*(xa[i] + xa_last[i]);
614: *changed_x = PETSC_TRUE;
615: 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]));
616: }
617: }
618: DMDAVecRestoreArray(da,check->last_step,&xa_last);
619: DMDAVecRestoreArray(da,x,&xa);
620: }
621: VecCopy(x,check->last_step);
622: return(0);
623: }
625: /* ------------------------------------------------------------------- */
626: /*
627: PostSetSubKSP - Optional user-defined routine that reset SubKSP options when hierarchical bjacobi PC is used
628: e.g,
629: 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
630: Set by SNESLineSearchSetPostCheck().
632: Input Parameters:
633: linesearch - the LineSearch context
634: xcurrent - current solution
635: y - search direction and length
636: x - the new candidate iterate
638: Output Parameters:
639: y - proposed step (search direction and length) (possibly changed)
640: x - current iterate (possibly modified)
642: */
643: PetscErrorCode PostSetSubKSP(SNESLineSearch linesearch,Vec xcurrent,Vec y,Vec x,PetscBool *changed_y,PetscBool *changed_x, void * ctx)
644: {
646: SetSubKSPCtx *check;
647: PetscInt iter,its,sub_its,maxit;
648: KSP ksp,sub_ksp,*sub_ksps;
649: PC pc;
650: PetscReal ksp_ratio;
651: SNES snes;
654: SNESLineSearchGetSNES(linesearch, &snes);
655: check = (SetSubKSPCtx*)ctx;
656: SNESGetIterationNumber(snes,&iter);
657: SNESGetKSP(snes,&ksp);
658: KSPGetPC(ksp,&pc);
659: PCBJacobiGetSubKSP(pc,NULL,NULL,&sub_ksps);
660: sub_ksp = sub_ksps[0];
661: KSPGetIterationNumber(ksp,&its); /* outer KSP iteration number */
662: KSPGetIterationNumber(sub_ksp,&sub_its); /* inner KSP iteration number */
664: if (iter) {
665: PetscPrintf(PETSC_COMM_WORLD," ...PostCheck snes iteration %D, ksp_it %D %D, subksp_it %D\n",iter,check->its0,its,sub_its);
666: ksp_ratio = ((PetscReal)(its))/check->its0;
667: maxit = (PetscInt)(ksp_ratio*sub_its + 0.5);
668: if (maxit < 2) maxit = 2;
669: KSPSetTolerances(sub_ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,maxit);
670: PetscPrintf(PETSC_COMM_WORLD," ...ksp_ratio %g, new maxit %D\n\n",(double)ksp_ratio,maxit);
671: }
672: check->its0 = its; /* save current outer KSP iteration number */
673: return(0);
674: }
676: /* ------------------------------------------------------------------- */
677: /*
678: MatrixFreePreconditioner - This routine demonstrates the use of a
679: user-provided preconditioner. This code implements just the null
680: preconditioner, which of course is not recommended for general use.
682: Input Parameters:
683: + pc - preconditioner
684: - x - input vector
686: Output Parameter:
687: . y - preconditioned vector
688: */
689: PetscErrorCode MatrixFreePreconditioner(PC pc,Vec x,Vec y)
690: {
692: VecCopy(x,y);
693: return 0;
694: }
696: /*TEST
698: test:
699: args: -nox -snes_monitor_cancel -snes_monitor_short -ksp_gmres_cgs_refinement_type refine_always
701: test:
702: suffix: 2
703: nsize: 3
704: args: -nox -pc_type asm -mat_type mpiaij -snes_monitor_cancel -snes_monitor_short -ksp_gmres_cgs_refinement_type refine_always
706: test:
707: suffix: 3
708: nsize: 2
709: args: -nox -snes_monitor_cancel -snes_monitor_short -ksp_gmres_cgs_refinement_type refine_always
711: test:
712: suffix: 4
713: args: -nox -pre_check_iterates -post_check_iterates
715: test:
716: suffix: 5
717: requires: double !complex !single
718: nsize: 2
719: args: -nox -snes_test_jacobian -snes_test_jacobian_view
721: test:
722: suffix: 6
723: requires: double !complex !single
724: nsize: 4
725: args: -test_jacobian_domain_error -snes_converged_reason -snes_check_jacobian_domain_error 1
727: test:
728: suffix: 7
729: requires: double !complex !single
730: nsize: 4
731: args: -test_jacobian_domain_error -snes_converged_reason -snes_type newtontr -snes_check_jacobian_domain_error 1
733: test:
734: suffix: 8
735: requires: double !complex !single
736: nsize: 4
737: args: -test_jacobian_domain_error -snes_converged_reason -snes_type vinewtonrsls -snes_check_jacobian_domain_error 1
739: test:
740: suffix: 9
741: requires: double !complex !single
742: nsize: 4
743: args: -test_jacobian_domain_error -snes_converged_reason -snes_type vinewtonssls -snes_check_jacobian_domain_error 1
745: test:
746: suffix: 10
747: requires: double !complex !single
748: nsize: 4
749: args: -test_jacobian_domain_error -snes_converged_reason -snes_type qn -snes_qn_scale_type jacobian -snes_check_jacobian_domain_error 1
751: test:
752: suffix: 11
753: requires: double !complex !single
754: nsize: 4
755: 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
757: test:
758: suffix: 12
759: args: -view_initial
760: filter: grep -v "type:"
762: TEST*/