Actual source code: itfunc.c
petsc-3.5.4 2015-05-23
2: /*
3: Interface KSP routines that the user calls.
4: */
6: #include <petsc-private/kspimpl.h> /*I "petscksp.h" I*/
7: #include <petscdm.h>
11: /*@
12: KSPComputeExtremeSingularValues - Computes the extreme singular values
13: for the preconditioned operator. Called after or during KSPSolve().
15: Not Collective
17: Input Parameter:
18: . ksp - iterative context obtained from KSPCreate()
20: Output Parameters:
21: . emin, emax - extreme singular values
23: Options Database Keys:
24: . -ksp_compute_singularvalues - compute extreme singular values and print when KSPSolve completes.
26: Notes:
27: One must call KSPSetComputeSingularValues() before calling KSPSetUp()
28: (or use the option -ksp_compute_eigenvalues) in order for this routine to work correctly.
30: Many users may just want to use the monitoring routine
31: KSPMonitorSingularValue() (which can be set with option -ksp_monitor_singular_value)
32: to print the extreme singular values at each iteration of the linear solve.
34: Estimates of the smallest singular value may be very inaccurate, especially if the Krylov method has not converged.
35: The largest singular value is usually accurate to within a few percent if the method has converged, but is still not
36: intended for eigenanalysis.
38: Disable restarts if using KSPGMRES, otherwise this estimate will only be using those iterations after the last
39: restart. See KSPGMRESSetRestart() for more details.
41: Level: advanced
43: .keywords: KSP, compute, extreme, singular, values
45: .seealso: KSPSetComputeSingularValues(), KSPMonitorSingularValue(), KSPComputeEigenvalues()
46: @*/
47: PetscErrorCode KSPComputeExtremeSingularValues(KSP ksp,PetscReal *emax,PetscReal *emin)
48: {
55: if (!ksp->calc_sings) SETERRQ(PetscObjectComm((PetscObject)ksp),4,"Singular values not requested before KSPSetUp()");
57: if (ksp->ops->computeextremesingularvalues) {
58: (*ksp->ops->computeextremesingularvalues)(ksp,emax,emin);
59: } else {
60: *emin = -1.0;
61: *emax = -1.0;
62: }
63: return(0);
64: }
68: /*@
69: KSPComputeEigenvalues - Computes the extreme eigenvalues for the
70: preconditioned operator. Called after or during KSPSolve().
72: Not Collective
74: Input Parameter:
75: + ksp - iterative context obtained from KSPCreate()
76: - n - size of arrays r and c. The number of eigenvalues computed (neig) will, in
77: general, be less than this.
79: Output Parameters:
80: + r - real part of computed eigenvalues, provided by user with a dimension of at least n
81: . c - complex part of computed eigenvalues, provided by user with a dimension of at least n
82: - neig - actual number of eigenvalues computed (will be less than or equal to n)
84: Options Database Keys:
85: + -ksp_compute_eigenvalues - Prints eigenvalues to stdout
86: - -ksp_plot_eigenvalues - Plots eigenvalues in an x-window display
88: Notes:
89: The number of eigenvalues estimated depends on the size of the Krylov space
90: generated during the KSPSolve() ; for example, with
91: CG it corresponds to the number of CG iterations, for GMRES it is the number
92: of GMRES iterations SINCE the last restart. Any extra space in r[] and c[]
93: will be ignored.
95: KSPComputeEigenvalues() does not usually provide accurate estimates; it is
96: intended only for assistance in understanding the convergence of iterative
97: methods, not for eigenanalysis. For accurate computation of eigenvalues we recommend using
98: the excellant package SLEPc.
100: One must call KSPSetComputeEigenvalues() before calling KSPSetUp()
101: in order for this routine to work correctly.
103: Many users may just want to use the monitoring routine
104: KSPMonitorSingularValue() (which can be set with option -ksp_monitor_singular_value)
105: to print the singular values at each iteration of the linear solve.
107: Level: advanced
109: .keywords: KSP, compute, extreme, singular, values
111: .seealso: KSPSetComputeSingularValues(), KSPMonitorSingularValue(), KSPComputeExtremeSingularValues()
112: @*/
113: PetscErrorCode KSPComputeEigenvalues(KSP ksp,PetscInt n,PetscReal r[],PetscReal c[],PetscInt *neig)
114: {
121: if (n<0) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_OUTOFRANGE,"Requested < 0 Eigenvalues");
123: if (!ksp->calc_sings) SETERRQ(PetscObjectComm((PetscObject)ksp),4,"Eigenvalues not requested before KSPSetUp()");
125: if (ksp->ops->computeeigenvalues) {
126: (*ksp->ops->computeeigenvalues)(ksp,n,r,c,neig);
127: } else {
128: *neig = 0;
129: }
130: return(0);
131: }
135: /*@
136: KSPSetUpOnBlocks - Sets up the preconditioner for each block in
137: the block Jacobi, block Gauss-Seidel, and overlapping Schwarz
138: methods.
140: Collective on KSP
142: Input Parameter:
143: . ksp - the KSP context
145: Notes:
146: KSPSetUpOnBlocks() is a routine that the user can optinally call for
147: more precise profiling (via -log_summary) of the setup phase for these
148: block preconditioners. If the user does not call KSPSetUpOnBlocks(),
149: it will automatically be called from within KSPSolve().
151: Calling KSPSetUpOnBlocks() is the same as calling PCSetUpOnBlocks()
152: on the PC context within the KSP context.
154: Level: advanced
156: .keywords: KSP, setup, blocks
158: .seealso: PCSetUpOnBlocks(), KSPSetUp(), PCSetUp()
159: @*/
160: PetscErrorCode KSPSetUpOnBlocks(KSP ksp)
161: {
166: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
167: PCSetUpOnBlocks(ksp->pc);
168: return(0);
169: }
173: /*@
174: KSPSetReusePreconditioner - reuse the current preconditioner, do not construct a new one even if the operator changes
176: Collective on KSP
178: Input Parameters:
179: + ksp - iterative context obtained from KSPCreate()
180: - flag - PETSC_TRUE to reuse the current preconditioner
182: Level: intermediate
184: .keywords: KSP, setup
186: .seealso: KSPCreate(), KSPSolve(), KSPDestroy(), PCSetReusePreconditioner()
187: @*/
188: PetscErrorCode KSPSetReusePreconditioner(KSP ksp,PetscBool flag)
189: {
194: PCSetReusePreconditioner(ksp->pc,flag);
195: return(0);
196: }
200: /*@
201: KSPSetUp - Sets up the internal data structures for the
202: later use of an iterative solver.
204: Collective on KSP
206: Input Parameter:
207: . ksp - iterative context obtained from KSPCreate()
209: Level: developer
211: .keywords: KSP, setup
213: .seealso: KSPCreate(), KSPSolve(), KSPDestroy()
214: @*/
215: PetscErrorCode KSPSetUp(KSP ksp)
216: {
218: Mat A,B;
223: /* reset the convergence flag from the previous solves */
224: ksp->reason = KSP_CONVERGED_ITERATING;
226: if (!((PetscObject)ksp)->type_name) {
227: KSPSetType(ksp,KSPGMRES);
228: }
229: KSPSetUpNorms_Private(ksp,&ksp->normtype,&ksp->pc_side);
231: if (ksp->dmActive && !ksp->setupstage) {
232: /* first time in so build matrix and vector data structures using DM */
233: if (!ksp->vec_rhs) {DMCreateGlobalVector(ksp->dm,&ksp->vec_rhs);}
234: if (!ksp->vec_sol) {DMCreateGlobalVector(ksp->dm,&ksp->vec_sol);}
235: DMCreateMatrix(ksp->dm,&A);
236: KSPSetOperators(ksp,A,A);
237: PetscObjectDereference((PetscObject)A);
238: }
240: if (ksp->dmActive) {
241: DMKSP kdm;
242: DMGetDMKSP(ksp->dm,&kdm);
244: if (kdm->ops->computeinitialguess && ksp->setupstage != KSP_SETUP_NEWRHS) {
245: /* only computes initial guess the first time through */
246: (*kdm->ops->computeinitialguess)(ksp,ksp->vec_sol,kdm->initialguessctx);
247: KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
248: }
249: if (kdm->ops->computerhs) {
250: (*kdm->ops->computerhs)(ksp,ksp->vec_rhs,kdm->rhsctx);
251: }
253: if (ksp->setupstage != KSP_SETUP_NEWRHS) {
254: if (kdm->ops->computeoperators) {
255: KSPGetOperators(ksp,&A,&B);
256: (*kdm->ops->computeoperators)(ksp,A,B,kdm->operatorsctx);
257: KSPSetOperators(ksp,A,B);
258: } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"You called KSPSetDM() but did not use DMKSPSetComputeOperators() or KSPSetDMActive(dm,PETSC_FALSE);");
259: }
260: }
262: if (ksp->setupstage == KSP_SETUP_NEWRHS) return(0);
263: PetscLogEventBegin(KSP_SetUp,ksp,ksp->vec_rhs,ksp->vec_sol,0);
265: switch (ksp->setupstage) {
266: case KSP_SETUP_NEW:
267: (*ksp->ops->setup)(ksp);
268: break;
269: case KSP_SETUP_NEWMATRIX: { /* This should be replaced with a more general mechanism */
270: KSPChebyshevSetNewMatrix(ksp);
271: } break;
272: default: break;
273: }
275: /* scale the matrix if requested */
276: if (ksp->dscale) {
277: Mat mat,pmat;
278: PetscScalar *xx;
279: PetscInt i,n;
280: PetscBool zeroflag = PETSC_FALSE;
281: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
282: PCGetOperators(ksp->pc,&mat,&pmat);
283: if (!ksp->diagonal) { /* allocate vector to hold diagonal */
284: MatGetVecs(pmat,&ksp->diagonal,0);
285: }
286: MatGetDiagonal(pmat,ksp->diagonal);
287: VecGetLocalSize(ksp->diagonal,&n);
288: VecGetArray(ksp->diagonal,&xx);
289: for (i=0; i<n; i++) {
290: if (xx[i] != 0.0) xx[i] = 1.0/PetscSqrtReal(PetscAbsScalar(xx[i]));
291: else {
292: xx[i] = 1.0;
293: zeroflag = PETSC_TRUE;
294: }
295: }
296: VecRestoreArray(ksp->diagonal,&xx);
297: if (zeroflag) {
298: PetscInfo(ksp,"Zero detected in diagonal of matrix, using 1 at those locations\n");
299: }
300: MatDiagonalScale(pmat,ksp->diagonal,ksp->diagonal);
301: if (mat != pmat) {MatDiagonalScale(mat,ksp->diagonal,ksp->diagonal);}
302: ksp->dscalefix2 = PETSC_FALSE;
303: }
304: PetscLogEventEnd(KSP_SetUp,ksp,ksp->vec_rhs,ksp->vec_sol,0);
305: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
306: PCSetUp(ksp->pc);
307: if (ksp->nullsp) {
308: PetscBool test = PETSC_FALSE;
309: PetscOptionsGetBool(((PetscObject)ksp)->prefix,"-ksp_test_null_space",&test,NULL);
310: if (test) {
311: Mat mat;
312: PCGetOperators(ksp->pc,&mat,NULL);
313: MatNullSpaceTest(ksp->nullsp,mat,NULL);
314: }
315: }
316: ksp->setupstage = KSP_SETUP_NEWRHS;
317: return(0);
318: }
320: #include <petscdraw.h>
323: /*@
324: KSPSolve - Solves linear system.
326: Collective on KSP
328: Parameter:
329: + ksp - iterative context obtained from KSPCreate()
330: . b - the right hand side vector
331: - x - the solution (this may be the same vector as b, then b will be overwritten with answer)
333: Options Database Keys:
334: + -ksp_compute_eigenvalues - compute preconditioned operators eigenvalues
335: . -ksp_plot_eigenvalues - plot the computed eigenvalues in an X-window
336: . -ksp_compute_eigenvalues_explicitly - compute the eigenvalues by forming the dense operator and useing LAPACK
337: . -ksp_plot_eigenvalues_explicitly - plot the explicitly computing eigenvalues
338: . -ksp_view_mat binary - save matrix to the default binary viewer
339: . -ksp_view_pmat binary - save matrix used to build preconditioner to the default binary viewer
340: . -ksp_view_rhs binary - save right hand side vector to the default binary viewer
341: . -ksp_view_solution binary - save computed solution vector to the default binary viewer
342: (can be read later with src/ksp/examples/tutorials/ex10.c for testing solvers)
343: . -ksp_view_mat_explicit - for matrix-free operators, computes the matrix entries and views them
344: . -ksp_view_preconditioned_operator_explicit - computes the product of the preconditioner and matrix as an explicit matrix and views it
345: . -ksp_converged_reason - print reason for converged or diverged, also prints number of iterations
346: . -ksp_final_residual - print 2-norm of true linear system residual at the end of the solution process
347: - -ksp_view - print the ksp data structure at the end of the system solution
349: Notes:
351: If one uses KSPSetDM() then x or b need not be passed. Use KSPGetSolution() to access the solution in this case.
353: The operator is specified with KSPSetOperators().
355: Call KSPGetConvergedReason() to determine if the solver converged or failed and
356: why. The number of iterations can be obtained from KSPGetIterationNumber().
358: If using a direct method (e.g., via the KSP solver
359: KSPPREONLY and a preconditioner such as PCLU/PCILU),
360: then its=1. See KSPSetTolerances() and KSPConvergedDefault()
361: for more details.
363: Understanding Convergence:
364: The routines KSPMonitorSet(), KSPComputeEigenvalues(), and
365: KSPComputeEigenvaluesExplicitly() provide information on additional
366: options to monitor convergence and print eigenvalue information.
368: Level: beginner
370: .keywords: KSP, solve, linear system
372: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPConvergedDefault(),
373: KSPSolveTranspose(), KSPGetIterationNumber()
374: @*/
375: PetscErrorCode KSPSolve(KSP ksp,Vec b,Vec x)
376: {
377: PetscErrorCode ierr;
378: PetscMPIInt rank;
379: PetscBool flag1,flag2,flag3,flg = PETSC_FALSE,inXisinB=PETSC_FALSE,guess_zero;
380: Mat mat,premat;
387: if (x && x == b) {
388: if (!ksp->guess_zero) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"Cannot use x == b with nonzero initial guess");
389: VecDuplicate(b,&x);
390: inXisinB = PETSC_TRUE;
391: }
392: if (b) {
393: PetscObjectReference((PetscObject)b);
394: VecDestroy(&ksp->vec_rhs);
395: ksp->vec_rhs = b;
396: }
397: if (x) {
398: PetscObjectReference((PetscObject)x);
399: VecDestroy(&ksp->vec_sol);
400: ksp->vec_sol = x;
401: }
402: KSPViewFromOptions(ksp,NULL,"-ksp_view_pre");
404: if (ksp->presolve) {
405: (*ksp->presolve)(ksp,ksp->vec_rhs,ksp->vec_sol,ksp->prectx);
406: }
407: PetscLogEventBegin(KSP_Solve,ksp,ksp->vec_rhs,ksp->vec_sol,0);
409: /* reset the residual history list if requested */
410: if (ksp->res_hist_reset) ksp->res_hist_len = 0;
411: ksp->transpose_solve = PETSC_FALSE;
413: if (ksp->guess) {
414: KSPFischerGuessFormGuess(ksp->guess,ksp->vec_rhs,ksp->vec_sol);
415: ksp->guess_zero = PETSC_FALSE;
416: }
417: /* KSPSetUp() scales the matrix if needed */
418: KSPSetUp(ksp);
419: KSPSetUpOnBlocks(ksp);
421: /* diagonal scale RHS if called for */
422: if (ksp->dscale) {
423: VecPointwiseMult(ksp->vec_rhs,ksp->vec_rhs,ksp->diagonal);
424: /* second time in, but matrix was scaled back to original */
425: if (ksp->dscalefix && ksp->dscalefix2) {
426: Mat mat,pmat;
428: PCGetOperators(ksp->pc,&mat,&pmat);
429: MatDiagonalScale(pmat,ksp->diagonal,ksp->diagonal);
430: if (mat != pmat) {MatDiagonalScale(mat,ksp->diagonal,ksp->diagonal);}
431: }
433: /* scale initial guess */
434: if (!ksp->guess_zero) {
435: if (!ksp->truediagonal) {
436: VecDuplicate(ksp->diagonal,&ksp->truediagonal);
437: VecCopy(ksp->diagonal,ksp->truediagonal);
438: VecReciprocal(ksp->truediagonal);
439: }
440: VecPointwiseMult(ksp->vec_sol,ksp->vec_sol,ksp->truediagonal);
441: }
442: }
443: PCPreSolve(ksp->pc,ksp);
445: if (ksp->guess_zero) { VecSet(ksp->vec_sol,0.0);}
446: if (ksp->guess_knoll) {
447: PCApply(ksp->pc,ksp->vec_rhs,ksp->vec_sol);
448: KSP_RemoveNullSpace(ksp,ksp->vec_sol);
449: ksp->guess_zero = PETSC_FALSE;
450: }
452: /* can we mark the initial guess as zero for this solve? */
453: guess_zero = ksp->guess_zero;
454: if (!ksp->guess_zero) {
455: PetscReal norm;
457: VecNormAvailable(ksp->vec_sol,NORM_2,&flg,&norm);
458: if (flg && !norm) ksp->guess_zero = PETSC_TRUE;
459: }
460: (*ksp->ops->solve)(ksp);
461: ksp->guess_zero = guess_zero;
463: if (!ksp->reason) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_PLIB,"Internal error, solver returned without setting converged reason");
464: if (ksp->printreason) {
465: PetscViewerASCIIAddTab(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)),((PetscObject)ksp)->tablevel);
466: if (ksp->reason > 0) {
467: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)),"Linear solve converged due to %s iterations %D\n",KSPConvergedReasons[ksp->reason],ksp->its);
468: } else {
469: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)),"Linear solve did not converge due to %s iterations %D\n",KSPConvergedReasons[ksp->reason],ksp->its);
470: }
471: PetscViewerASCIISubtractTab(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)ksp)),((PetscObject)ksp)->tablevel);
472: }
473: PCPostSolve(ksp->pc,ksp);
475: /* diagonal scale solution if called for */
476: if (ksp->dscale) {
477: VecPointwiseMult(ksp->vec_sol,ksp->vec_sol,ksp->diagonal);
478: /* unscale right hand side and matrix */
479: if (ksp->dscalefix) {
480: Mat mat,pmat;
482: VecReciprocal(ksp->diagonal);
483: VecPointwiseMult(ksp->vec_rhs,ksp->vec_rhs,ksp->diagonal);
484: PCGetOperators(ksp->pc,&mat,&pmat);
485: MatDiagonalScale(pmat,ksp->diagonal,ksp->diagonal);
486: if (mat != pmat) {MatDiagonalScale(mat,ksp->diagonal,ksp->diagonal);}
487: VecReciprocal(ksp->diagonal);
488: ksp->dscalefix2 = PETSC_TRUE;
489: }
490: }
491: PetscLogEventEnd(KSP_Solve,ksp,ksp->vec_rhs,ksp->vec_sol,0);
492: if (ksp->postsolve) {
493: (*ksp->postsolve)(ksp,ksp->vec_rhs,ksp->vec_sol,ksp->postctx);
494: }
496: if (ksp->guess) {
497: KSPFischerGuessUpdate(ksp->guess,ksp->vec_sol);
498: }
500: MPI_Comm_rank(PetscObjectComm((PetscObject)ksp),&rank);
502: PCGetOperators(ksp->pc,&mat,&premat);
503: MatViewFromOptions(mat,((PetscObject)ksp)->prefix,"-ksp_view_mat");
504: MatViewFromOptions(premat,((PetscObject)ksp)->prefix,"-ksp_view_pmat");
505: VecViewFromOptions(ksp->vec_rhs,((PetscObject)ksp)->prefix,"-ksp_view_rhs");
507: flag1 = PETSC_FALSE;
508: flag2 = PETSC_FALSE;
509: flag3 = PETSC_FALSE;
510: PetscOptionsGetBool(((PetscObject)ksp)->prefix,"-ksp_compute_eigenvalues",&flag1,NULL);
511: PetscOptionsGetBool(((PetscObject)ksp)->prefix,"-ksp_plot_eigenvalues",&flag2,NULL);
512: PetscOptionsGetBool(((PetscObject)ksp)->prefix,"-ksp_plot_eigencontours",&flag3,NULL);
513: if (flag1 || flag2 || flag3) {
514: PetscInt nits,n,i,neig;
515: PetscReal *r,*c;
517: KSPGetIterationNumber(ksp,&nits);
518: n = nits+2;
520: if (!nits) {
521: PetscPrintf(PetscObjectComm((PetscObject)ksp),"Zero iterations in solver, cannot approximate any eigenvalues\n");
522: } else {
523: PetscMalloc2(n,&r,n,&c);
524: KSPComputeEigenvalues(ksp,n,r,c,&neig);
525: if (flag1) {
526: PetscPrintf(PetscObjectComm((PetscObject)ksp),"Iteratively computed eigenvalues\n");
527: for (i=0; i<neig; i++) {
528: if (c[i] >= 0.0) {
529: PetscPrintf(PetscObjectComm((PetscObject)ksp),"%g + %gi\n",(double)r[i],(double)c[i]);
530: } else {
531: PetscPrintf(PetscObjectComm((PetscObject)ksp),"%g - %gi\n",(double)r[i],-(double)c[i]);
532: }
533: }
534: }
535: if (flag2 && !rank) {
536: PetscDraw draw;
537: PetscDrawSP drawsp;
539: if (!ksp->eigviewer) {
540: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Iteratively Computed Eigenvalues",PETSC_DECIDE,PETSC_DECIDE,400,400,&ksp->eigviewer);
541: }
542: PetscViewerDrawGetDraw(ksp->eigviewer,0,&draw);
543: PetscDrawSPCreate(draw,1,&drawsp);
544: PetscDrawSPReset(drawsp);
545: for (i=0; i<neig; i++) {
546: PetscDrawSPAddPoint(drawsp,r+i,c+i);
547: }
548: PetscDrawSPDraw(drawsp,PETSC_TRUE);
549: PetscDrawSPDestroy(&drawsp);
550: }
551: if (flag3 && !rank) {
552: KSPPlotEigenContours_Private(ksp,neig,r,c);
553: }
554: PetscFree2(r,c);
555: }
556: }
558: flag1 = PETSC_FALSE;
559: PetscOptionsGetBool(((PetscObject)ksp)->prefix,"-ksp_compute_singularvalues",&flag1,NULL);
560: if (flag1) {
561: PetscInt nits;
563: KSPGetIterationNumber(ksp,&nits);
564: if (!nits) {
565: PetscPrintf(PetscObjectComm((PetscObject)ksp),"Zero iterations in solver, cannot approximate any singular values\n");
566: } else {
567: PetscReal emax,emin;
569: KSPComputeExtremeSingularValues(ksp,&emax,&emin);
570: PetscPrintf(PetscObjectComm((PetscObject)ksp),"Iteratively computed extreme singular values: max %g min %g max/min %g\n",(double)emax,(double)emin,(double)(emax/emin));
571: }
572: }
575: flag1 = PETSC_FALSE;
576: flag2 = PETSC_FALSE;
577: PetscOptionsGetBool(((PetscObject)ksp)->prefix,"-ksp_compute_eigenvalues_explicitly",&flag1,NULL);
578: PetscOptionsGetBool(((PetscObject)ksp)->prefix,"-ksp_plot_eigenvalues_explicitly",&flag2,NULL);
579: if (flag1 || flag2) {
580: PetscInt n,i;
581: PetscReal *r,*c;
582: VecGetSize(ksp->vec_sol,&n);
583: PetscMalloc2(n,&r,n,&c);
584: KSPComputeEigenvaluesExplicitly(ksp,n,r,c);
585: if (flag1) {
586: PetscPrintf(PetscObjectComm((PetscObject)ksp),"Explicitly computed eigenvalues\n");
587: for (i=0; i<n; i++) {
588: if (c[i] >= 0.0) {
589: PetscPrintf(PetscObjectComm((PetscObject)ksp),"%g + %gi\n",(double)r[i],(double)c[i]);
590: } else {
591: PetscPrintf(PetscObjectComm((PetscObject)ksp),"%g - %gi\n",(double)r[i],-(double)c[i]);
592: }
593: }
594: }
595: if (flag2 && !rank) {
596: PetscDraw draw;
597: PetscDrawSP drawsp;
599: if (!ksp->eigviewer) {
600: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Explicitly Computed Eigenvalues",0,320,400,400,&ksp->eigviewer);
601: }
602: PetscViewerDrawGetDraw(ksp->eigviewer,0,&draw);
603: PetscDrawSPCreate(draw,1,&drawsp);
604: PetscDrawSPReset(drawsp);
605: for (i=0; i<n; i++) {
606: PetscDrawSPAddPoint(drawsp,r+i,c+i);
607: }
608: PetscDrawSPDraw(drawsp,PETSC_TRUE);
609: PetscDrawSPDestroy(&drawsp);
610: }
611: PetscFree2(r,c);
612: }
614: PetscOptionsHasName(((PetscObject)ksp)->prefix,"-ksp_view_mat_explicit",&flag2);
615: if (flag2) {
616: Mat A,B;
617: PCGetOperators(ksp->pc,&A,NULL);
618: MatComputeExplicitOperator(A,&B);
619: MatViewFromOptions(B,((PetscObject)ksp)->prefix,"-ksp_view_mat_explicit");
620: MatDestroy(&B);
621: }
622: PetscOptionsHasName(((PetscObject)ksp)->prefix,"-ksp_view_preconditioned_operator_explicit",&flag2);
623: if (flag2) {
624: Mat B;
625: KSPComputeExplicitOperator(ksp,&B);
626: MatViewFromOptions(B,((PetscObject)ksp)->prefix,"-ksp_view_preconditioned_operator_explicit");
627: MatDestroy(&B);
628: }
629: KSPViewFromOptions(ksp,NULL,"-ksp_view");
631: flg = PETSC_FALSE;
632: PetscOptionsGetBool(((PetscObject)ksp)->prefix,"-ksp_final_residual",&flg,NULL);
633: if (flg) {
634: Mat A;
635: Vec t;
636: PetscReal norm;
637: if (ksp->dscale && !ksp->dscalefix) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"Cannot compute final scale with -ksp_diagonal_scale except also with -ksp_diagonal_scale_fix");
638: PCGetOperators(ksp->pc,&A,NULL);
639: VecDuplicate(ksp->vec_rhs,&t);
640: KSP_MatMult(ksp,A,ksp->vec_sol,t);
641: VecAYPX(t, -1.0, ksp->vec_rhs);
642: VecNorm(t,NORM_2,&norm);
643: VecDestroy(&t);
644: PetscPrintf(PetscObjectComm((PetscObject)ksp),"KSP final norm of residual %g\n",(double)norm);
645: }
646: VecViewFromOptions(ksp->vec_sol,((PetscObject)ksp)->prefix,"-ksp_view_solution");
648: if (inXisinB) {
649: VecCopy(x,b);
650: VecDestroy(&x);
651: }
652: PetscObjectSAWsBlock((PetscObject)ksp);
653: if (ksp->errorifnotconverged && ksp->reason < 0) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_NOT_CONVERGED,"KSPSolve has not converged");
654: return(0);
655: }
659: /*@
660: KSPSolveTranspose - Solves the transpose of a linear system.
662: Collective on KSP
664: Input Parameter:
665: + ksp - iterative context obtained from KSPCreate()
666: . b - right hand side vector
667: - x - solution vector
669: Notes: For complex numbers this solve the non-Hermitian transpose system.
671: Developer Notes: We need to implement a KSPSolveHermitianTranspose()
673: Level: developer
675: .keywords: KSP, solve, linear system
677: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPConvergedDefault(),
678: KSPSolve()
679: @*/
681: PetscErrorCode KSPSolveTranspose(KSP ksp,Vec b,Vec x)
682: {
684: PetscBool inXisinB=PETSC_FALSE;
690: if (x == b) {
691: VecDuplicate(b,&x);
692: inXisinB = PETSC_TRUE;
693: }
694: PetscObjectReference((PetscObject)b);
695: PetscObjectReference((PetscObject)x);
696: VecDestroy(&ksp->vec_rhs);
697: VecDestroy(&ksp->vec_sol);
699: ksp->vec_rhs = b;
700: ksp->vec_sol = x;
701: ksp->transpose_solve = PETSC_TRUE;
703: KSPSetUp(ksp);
704: if (ksp->guess_zero) { VecSet(ksp->vec_sol,0.0);}
705: (*ksp->ops->solve)(ksp);
706: if (!ksp->reason) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_PLIB,"Internal error, solver returned without setting converged reason");
707: if (ksp->printreason) {
708: if (ksp->reason > 0) {
709: PetscPrintf(PetscObjectComm((PetscObject)ksp),"Linear solve converged due to %s iterations %D\n",KSPConvergedReasons[ksp->reason],ksp->its);
710: } else {
711: PetscPrintf(PetscObjectComm((PetscObject)ksp),"Linear solve did not converge due to %s iterations %D\n",KSPConvergedReasons[ksp->reason],ksp->its);
712: }
713: }
714: if (inXisinB) {
715: VecCopy(x,b);
716: VecDestroy(&x);
717: }
718: if (ksp->errorifnotconverged && ksp->reason < 0) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_NOT_CONVERGED,"KSPSolve has not converged");
719: return(0);
720: }
724: /*@
725: KSPReset - Resets a KSP context to the kspsetupcalled = 0 state and removes any allocated Vecs and Mats
727: Collective on KSP
729: Input Parameter:
730: . ksp - iterative context obtained from KSPCreate()
732: Level: beginner
734: .keywords: KSP, destroy
736: .seealso: KSPCreate(), KSPSetUp(), KSPSolve()
737: @*/
738: PetscErrorCode KSPReset(KSP ksp)
739: {
744: if (!ksp) return(0);
745: if (ksp->ops->reset) {
746: (*ksp->ops->reset)(ksp);
747: }
748: if (ksp->pc) {PCReset(ksp->pc);}
749: KSPFischerGuessDestroy(&ksp->guess);
750: VecDestroyVecs(ksp->nwork,&ksp->work);
751: VecDestroy(&ksp->vec_rhs);
752: VecDestroy(&ksp->vec_sol);
753: VecDestroy(&ksp->diagonal);
754: VecDestroy(&ksp->truediagonal);
755: MatNullSpaceDestroy(&ksp->nullsp);
757: ksp->setupstage = KSP_SETUP_NEW;
758: return(0);
759: }
763: /*@
764: KSPDestroy - Destroys KSP context.
766: Collective on KSP
768: Input Parameter:
769: . ksp - iterative context obtained from KSPCreate()
771: Level: beginner
773: .keywords: KSP, destroy
775: .seealso: KSPCreate(), KSPSetUp(), KSPSolve()
776: @*/
777: PetscErrorCode KSPDestroy(KSP *ksp)
778: {
780: PC pc;
783: if (!*ksp) return(0);
785: if (--((PetscObject)(*ksp))->refct > 0) {*ksp = 0; return(0);}
787: PetscObjectSAWsViewOff((PetscObject)*ksp);
788: /*
789: Avoid a cascading call to PCReset(ksp->pc) from the following call:
790: PCReset() shouldn't be called from KSPDestroy() as it is unprotected by pc's
791: refcount (and may be shared, e.g., by other ksps).
792: */
793: pc = (*ksp)->pc;
794: (*ksp)->pc = NULL;
795: KSPReset((*ksp));
796: (*ksp)->pc = pc;
797: if ((*ksp)->ops->destroy) {(*(*ksp)->ops->destroy)(*ksp);}
799: DMDestroy(&(*ksp)->dm);
800: PCDestroy(&(*ksp)->pc);
801: PetscFree((*ksp)->res_hist_alloc);
802: if ((*ksp)->convergeddestroy) {
803: (*(*ksp)->convergeddestroy)((*ksp)->cnvP);
804: }
805: KSPMonitorCancel((*ksp));
806: PetscViewerDestroy(&(*ksp)->eigviewer);
807: PetscHeaderDestroy(ksp);
808: return(0);
809: }
813: /*@
814: KSPSetPCSide - Sets the preconditioning side.
816: Logically Collective on KSP
818: Input Parameter:
819: . ksp - iterative context obtained from KSPCreate()
821: Output Parameter:
822: . side - the preconditioning side, where side is one of
823: .vb
824: PC_LEFT - left preconditioning (default)
825: PC_RIGHT - right preconditioning
826: PC_SYMMETRIC - symmetric preconditioning
827: .ve
829: Options Database Keys:
830: . -ksp_pc_side <right,left,symmetric>
832: Notes:
833: Left preconditioning is used by default for most Krylov methods except KSPFGMRES which only supports right preconditioning.
834: Symmetric preconditioning is currently available only for the KSPQCG method. Note, however, that
835: symmetric preconditioning can be emulated by using either right or left
836: preconditioning and a pre or post processing step.
838: Setting the PC side often affects the default norm type. See KSPSetNormType() for details.
840: Level: intermediate
842: .keywords: KSP, set, right, left, symmetric, side, preconditioner, flag
844: .seealso: KSPGetPCSide(), KSPSetNormType()
845: @*/
846: PetscErrorCode KSPSetPCSide(KSP ksp,PCSide side)
847: {
851: ksp->pc_side = ksp->pc_side_set = side;
852: return(0);
853: }
857: /*@
858: KSPGetPCSide - Gets the preconditioning side.
860: Not Collective
862: Input Parameter:
863: . ksp - iterative context obtained from KSPCreate()
865: Output Parameter:
866: . side - the preconditioning side, where side is one of
867: .vb
868: PC_LEFT - left preconditioning (default)
869: PC_RIGHT - right preconditioning
870: PC_SYMMETRIC - symmetric preconditioning
871: .ve
873: Level: intermediate
875: .keywords: KSP, get, right, left, symmetric, side, preconditioner, flag
877: .seealso: KSPSetPCSide()
878: @*/
879: PetscErrorCode KSPGetPCSide(KSP ksp,PCSide *side)
880: {
886: KSPSetUpNorms_Private(ksp,&ksp->normtype,&ksp->pc_side);
887: *side = ksp->pc_side;
888: return(0);
889: }
893: /*@
894: KSPGetTolerances - Gets the relative, absolute, divergence, and maximum
895: iteration tolerances used by the default KSP convergence tests.
897: Not Collective
899: Input Parameter:
900: . ksp - the Krylov subspace context
902: Output Parameters:
903: + rtol - the relative convergence tolerance
904: . abstol - the absolute convergence tolerance
905: . dtol - the divergence tolerance
906: - maxits - maximum number of iterations
908: Notes:
909: The user can specify NULL for any parameter that is not needed.
911: Level: intermediate
913: .keywords: KSP, get, tolerance, absolute, relative, divergence, convergence,
914: maximum, iterations
916: .seealso: KSPSetTolerances()
917: @*/
918: PetscErrorCode KSPGetTolerances(KSP ksp,PetscReal *rtol,PetscReal *abstol,PetscReal *dtol,PetscInt *maxits)
919: {
922: if (abstol) *abstol = ksp->abstol;
923: if (rtol) *rtol = ksp->rtol;
924: if (dtol) *dtol = ksp->divtol;
925: if (maxits) *maxits = ksp->max_it;
926: return(0);
927: }
931: /*@
932: KSPSetTolerances - Sets the relative, absolute, divergence, and maximum
933: iteration tolerances used by the default KSP convergence testers.
935: Logically Collective on KSP
937: Input Parameters:
938: + ksp - the Krylov subspace context
939: . rtol - the relative convergence tolerance, relative decrease in the (possibly preconditioned) residual norm
940: . abstol - the absolute convergence tolerance absolute size of the (possibly preconditioned) residual norm
941: . dtol - the divergence tolerance, amount (possibly preconditioned) residual norm can increase before KSPConvergedDefault() concludes that the method is diverging
942: - maxits - maximum number of iterations to use
944: Options Database Keys:
945: + -ksp_atol <abstol> - Sets abstol
946: . -ksp_rtol <rtol> - Sets rtol
947: . -ksp_divtol <dtol> - Sets dtol
948: - -ksp_max_it <maxits> - Sets maxits
950: Notes:
951: Use PETSC_DEFAULT to retain the default value of any of the tolerances.
953: See KSPConvergedDefault() for details how these parameters are used in the default convergence test. See also KSPSetConvergenceTest()
954: for setting user-defined stopping criteria.
956: Level: intermediate
958: .keywords: KSP, set, tolerance, absolute, relative, divergence,
959: convergence, maximum, iterations
961: .seealso: KSPGetTolerances(), KSPConvergedDefault(), KSPSetConvergenceTest()
962: @*/
963: PetscErrorCode KSPSetTolerances(KSP ksp,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt maxits)
964: {
972: if (rtol != PETSC_DEFAULT) {
973: if (rtol < 0.0 || 1.0 <= rtol) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_OUTOFRANGE,"Relative tolerance %g must be non-negative and less than 1.0",(double)rtol);
974: ksp->rtol = rtol;
975: }
976: if (abstol != PETSC_DEFAULT) {
977: if (abstol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %g must be non-negative",(double)abstol);
978: ksp->abstol = abstol;
979: }
980: if (dtol != PETSC_DEFAULT) {
981: if (dtol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_OUTOFRANGE,"Divergence tolerance %g must be larger than 1.0",(double)dtol);
982: ksp->divtol = dtol;
983: }
984: if (maxits != PETSC_DEFAULT) {
985: if (maxits < 0) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",maxits);
986: ksp->max_it = maxits;
987: }
988: return(0);
989: }
993: /*@
994: KSPSetInitialGuessNonzero - Tells the iterative solver that the
995: initial guess is nonzero; otherwise KSP assumes the initial guess
996: is to be zero (and thus zeros it out before solving).
998: Logically Collective on KSP
1000: Input Parameters:
1001: + ksp - iterative context obtained from KSPCreate()
1002: - flg - PETSC_TRUE indicates the guess is non-zero, PETSC_FALSE indicates the guess is zero
1004: Options database keys:
1005: . -ksp_initial_guess_nonzero : use nonzero initial guess; this takes an optional truth value (0/1/no/yes/true/false)
1007: Level: beginner
1009: Notes:
1010: If this is not called the X vector is zeroed in the call to KSPSolve().
1012: .keywords: KSP, set, initial guess, nonzero
1014: .seealso: KSPGetInitialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
1015: @*/
1016: PetscErrorCode KSPSetInitialGuessNonzero(KSP ksp,PetscBool flg)
1017: {
1021: ksp->guess_zero = (PetscBool) !(int)flg;
1022: return(0);
1023: }
1027: /*@
1028: KSPGetInitialGuessNonzero - Determines whether the KSP solver is using
1029: a zero initial guess.
1031: Not Collective
1033: Input Parameter:
1034: . ksp - iterative context obtained from KSPCreate()
1036: Output Parameter:
1037: . flag - PETSC_TRUE if guess is nonzero, else PETSC_FALSE
1039: Level: intermediate
1041: .keywords: KSP, set, initial guess, nonzero
1043: .seealso: KSPSetInitialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
1044: @*/
1045: PetscErrorCode KSPGetInitialGuessNonzero(KSP ksp,PetscBool *flag)
1046: {
1050: if (ksp->guess_zero) *flag = PETSC_FALSE;
1051: else *flag = PETSC_TRUE;
1052: return(0);
1053: }
1057: /*@
1058: KSPSetErrorIfNotConverged - Causes KSPSolve() to generate an error if the solver has not converged.
1060: Logically Collective on KSP
1062: Input Parameters:
1063: + ksp - iterative context obtained from KSPCreate()
1064: - flg - PETSC_TRUE indicates you want the error generated
1066: Options database keys:
1067: . -ksp_error_if_not_converged : this takes an optional truth value (0/1/no/yes/true/false)
1069: Level: intermediate
1071: Notes:
1072: Normally PETSc continues if a linear solver fails to converge, you can call KSPGetConvergedReason() after a KSPSolve()
1073: to determine if it has converged.
1075: .keywords: KSP, set, initial guess, nonzero
1077: .seealso: KSPGetInitialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll(), KSPGetErrorIfNotConverged()
1078: @*/
1079: PetscErrorCode KSPSetErrorIfNotConverged(KSP ksp,PetscBool flg)
1080: {
1084: ksp->errorifnotconverged = flg;
1085: return(0);
1086: }
1090: /*@
1091: KSPGetErrorIfNotConverged - Will KSPSolve() generate an error if the solver does not converge?
1093: Not Collective
1095: Input Parameter:
1096: . ksp - iterative context obtained from KSPCreate()
1098: Output Parameter:
1099: . flag - PETSC_TRUE if it will generate an error, else PETSC_FALSE
1101: Level: intermediate
1103: .keywords: KSP, set, initial guess, nonzero
1105: .seealso: KSPSetInitialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll(), KSPSetErrorIfNotConverged()
1106: @*/
1107: PetscErrorCode KSPGetErrorIfNotConverged(KSP ksp,PetscBool *flag)
1108: {
1112: *flag = ksp->errorifnotconverged;
1113: return(0);
1114: }
1118: /*@
1119: KSPSetInitialGuessKnoll - Tells the iterative solver to use PCApply(pc,b,..) to compute the initial guess (The Knoll trick)
1121: Logically Collective on KSP
1123: Input Parameters:
1124: + ksp - iterative context obtained from KSPCreate()
1125: - flg - PETSC_TRUE or PETSC_FALSE
1127: Level: advanced
1130: .keywords: KSP, set, initial guess, nonzero
1132: .seealso: KSPGetInitialGuessKnoll(), KSPSetInitialGuessNonzero(), KSPGetInitialGuessNonzero()
1133: @*/
1134: PetscErrorCode KSPSetInitialGuessKnoll(KSP ksp,PetscBool flg)
1135: {
1139: ksp->guess_knoll = flg;
1140: return(0);
1141: }
1145: /*@
1146: KSPGetInitialGuessKnoll - Determines whether the KSP solver is using the Knoll trick (using PCApply(pc,b,...) to compute
1147: the initial guess
1149: Not Collective
1151: Input Parameter:
1152: . ksp - iterative context obtained from KSPCreate()
1154: Output Parameter:
1155: . flag - PETSC_TRUE if using Knoll trick, else PETSC_FALSE
1157: Level: advanced
1159: .keywords: KSP, set, initial guess, nonzero
1161: .seealso: KSPSetInitialGuessKnoll(), KSPSetInitialGuessNonzero(), KSPGetInitialGuessNonzero()
1162: @*/
1163: PetscErrorCode KSPGetInitialGuessKnoll(KSP ksp,PetscBool *flag)
1164: {
1168: *flag = ksp->guess_knoll;
1169: return(0);
1170: }
1174: /*@
1175: KSPGetComputeSingularValues - Gets the flag indicating whether the extreme singular
1176: values will be calculated via a Lanczos or Arnoldi process as the linear
1177: system is solved.
1179: Not Collective
1181: Input Parameter:
1182: . ksp - iterative context obtained from KSPCreate()
1184: Output Parameter:
1185: . flg - PETSC_TRUE or PETSC_FALSE
1187: Options Database Key:
1188: . -ksp_monitor_singular_value - Activates KSPSetComputeSingularValues()
1190: Notes:
1191: Currently this option is not valid for all iterative methods.
1193: Many users may just want to use the monitoring routine
1194: KSPMonitorSingularValue() (which can be set with option -ksp_monitor_singular_value)
1195: to print the singular values at each iteration of the linear solve.
1197: Level: advanced
1199: .keywords: KSP, set, compute, singular values
1201: .seealso: KSPComputeExtremeSingularValues(), KSPMonitorSingularValue()
1202: @*/
1203: PetscErrorCode KSPGetComputeSingularValues(KSP ksp,PetscBool *flg)
1204: {
1208: *flg = ksp->calc_sings;
1209: return(0);
1210: }
1214: /*@
1215: KSPSetComputeSingularValues - Sets a flag so that the extreme singular
1216: values will be calculated via a Lanczos or Arnoldi process as the linear
1217: system is solved.
1219: Logically Collective on KSP
1221: Input Parameters:
1222: + ksp - iterative context obtained from KSPCreate()
1223: - flg - PETSC_TRUE or PETSC_FALSE
1225: Options Database Key:
1226: . -ksp_monitor_singular_value - Activates KSPSetComputeSingularValues()
1228: Notes:
1229: Currently this option is not valid for all iterative methods.
1231: Many users may just want to use the monitoring routine
1232: KSPMonitorSingularValue() (which can be set with option -ksp_monitor_singular_value)
1233: to print the singular values at each iteration of the linear solve.
1235: Level: advanced
1237: .keywords: KSP, set, compute, singular values
1239: .seealso: KSPComputeExtremeSingularValues(), KSPMonitorSingularValue()
1240: @*/
1241: PetscErrorCode KSPSetComputeSingularValues(KSP ksp,PetscBool flg)
1242: {
1246: ksp->calc_sings = flg;
1247: return(0);
1248: }
1252: /*@
1253: KSPGetComputeEigenvalues - Gets the flag indicating that the extreme eigenvalues
1254: values will be calculated via a Lanczos or Arnoldi process as the linear
1255: system is solved.
1257: Not Collective
1259: Input Parameter:
1260: . ksp - iterative context obtained from KSPCreate()
1262: Output Parameter:
1263: . flg - PETSC_TRUE or PETSC_FALSE
1265: Notes:
1266: Currently this option is not valid for all iterative methods.
1268: Level: advanced
1270: .keywords: KSP, set, compute, eigenvalues
1272: .seealso: KSPComputeEigenvalues(), KSPComputeEigenvaluesExplicitly()
1273: @*/
1274: PetscErrorCode KSPGetComputeEigenvalues(KSP ksp,PetscBool *flg)
1275: {
1279: *flg = ksp->calc_sings;
1280: return(0);
1281: }
1285: /*@
1286: KSPSetComputeEigenvalues - Sets a flag so that the extreme eigenvalues
1287: values will be calculated via a Lanczos or Arnoldi process as the linear
1288: system is solved.
1290: Logically Collective on KSP
1292: Input Parameters:
1293: + ksp - iterative context obtained from KSPCreate()
1294: - flg - PETSC_TRUE or PETSC_FALSE
1296: Notes:
1297: Currently this option is not valid for all iterative methods.
1299: Level: advanced
1301: .keywords: KSP, set, compute, eigenvalues
1303: .seealso: KSPComputeEigenvalues(), KSPComputeEigenvaluesExplicitly()
1304: @*/
1305: PetscErrorCode KSPSetComputeEigenvalues(KSP ksp,PetscBool flg)
1306: {
1310: ksp->calc_sings = flg;
1311: return(0);
1312: }
1316: /*@
1317: KSPGetRhs - Gets the right-hand-side vector for the linear system to
1318: be solved.
1320: Not Collective
1322: Input Parameter:
1323: . ksp - iterative context obtained from KSPCreate()
1325: Output Parameter:
1326: . r - right-hand-side vector
1328: Level: developer
1330: .keywords: KSP, get, right-hand-side, rhs
1332: .seealso: KSPGetSolution(), KSPSolve()
1333: @*/
1334: PetscErrorCode KSPGetRhs(KSP ksp,Vec *r)
1335: {
1339: *r = ksp->vec_rhs;
1340: return(0);
1341: }
1345: /*@
1346: KSPGetSolution - Gets the location of the solution for the
1347: linear system to be solved. Note that this may not be where the solution
1348: is stored during the iterative process; see KSPBuildSolution().
1350: Not Collective
1352: Input Parameters:
1353: . ksp - iterative context obtained from KSPCreate()
1355: Output Parameters:
1356: . v - solution vector
1358: Level: developer
1360: .keywords: KSP, get, solution
1362: .seealso: KSPGetRhs(), KSPBuildSolution(), KSPSolve()
1363: @*/
1364: PetscErrorCode KSPGetSolution(KSP ksp,Vec *v)
1365: {
1369: *v = ksp->vec_sol;
1370: return(0);
1371: }
1375: /*@
1376: KSPSetPC - Sets the preconditioner to be used to calculate the
1377: application of the preconditioner on a vector.
1379: Collective on KSP
1381: Input Parameters:
1382: + ksp - iterative context obtained from KSPCreate()
1383: - pc - the preconditioner object
1385: Notes:
1386: Use KSPGetPC() to retrieve the preconditioner context (for example,
1387: to free it at the end of the computations).
1389: Level: developer
1391: .keywords: KSP, set, precondition, Binv
1393: .seealso: KSPGetPC()
1394: @*/
1395: PetscErrorCode KSPSetPC(KSP ksp,PC pc)
1396: {
1403: PetscObjectReference((PetscObject)pc);
1404: PCDestroy(&ksp->pc);
1405: ksp->pc = pc;
1406: PetscLogObjectParent((PetscObject)ksp,(PetscObject)ksp->pc);
1407: return(0);
1408: }
1412: /*@
1413: KSPGetPC - Returns a pointer to the preconditioner context
1414: set with KSPSetPC().
1416: Not Collective
1418: Input Parameters:
1419: . ksp - iterative context obtained from KSPCreate()
1421: Output Parameter:
1422: . pc - preconditioner context
1424: Level: developer
1426: .keywords: KSP, get, preconditioner, Binv
1428: .seealso: KSPSetPC()
1429: @*/
1430: PetscErrorCode KSPGetPC(KSP ksp,PC *pc)
1431: {
1437: if (!ksp->pc) {
1438: PCCreate(PetscObjectComm((PetscObject)ksp),&ksp->pc);
1439: PetscObjectIncrementTabLevel((PetscObject)ksp->pc,(PetscObject)ksp,0);
1440: PetscLogObjectParent((PetscObject)ksp,(PetscObject)ksp->pc);
1441: }
1442: *pc = ksp->pc;
1443: return(0);
1444: }
1448: /*@
1449: KSPMonitor - runs the user provided monitor routines, if they exist
1451: Collective on KSP
1453: Input Parameters:
1454: + ksp - iterative context obtained from KSPCreate()
1455: . it - iteration number
1456: - rnorm - relative norm of the residual
1458: Notes:
1459: This routine is called by the KSP implementations.
1460: It does not typically need to be called by the user.
1462: Level: developer
1464: .seealso: KSPMonitorSet()
1465: @*/
1466: PetscErrorCode KSPMonitor(KSP ksp,PetscInt it,PetscReal rnorm)
1467: {
1468: PetscInt i, n = ksp->numbermonitors;
1472: for (i=0; i<n; i++) {
1473: (*ksp->monitor[i])(ksp,it,rnorm,ksp->monitorcontext[i]);
1474: }
1475: return(0);
1476: }
1480: /*@C
1481: KSPMonitorSet - Sets an ADDITIONAL function to be called at every iteration to monitor
1482: the residual/error etc.
1484: Logically Collective on KSP
1486: Input Parameters:
1487: + ksp - iterative context obtained from KSPCreate()
1488: . monitor - pointer to function (if this is NULL, it turns off monitoring
1489: . mctx - [optional] context for private data for the
1490: monitor routine (use NULL if no context is desired)
1491: - monitordestroy - [optional] routine that frees monitor context
1492: (may be NULL)
1494: Calling Sequence of monitor:
1495: $ monitor (KSP ksp, int it, PetscReal rnorm, void *mctx)
1497: + ksp - iterative context obtained from KSPCreate()
1498: . it - iteration number
1499: . rnorm - (estimated) 2-norm of (preconditioned) residual
1500: - mctx - optional monitoring context, as set by KSPMonitorSet()
1502: Options Database Keys:
1503: + -ksp_monitor - sets KSPMonitorDefault()
1504: . -ksp_monitor_true_residual - sets KSPMonitorTrueResidualNorm()
1505: . -ksp_monitor_max - sets KSPMonitorTrueResidualMaxNorm()
1506: . -ksp_monitor_lg_residualnorm - sets line graph monitor,
1507: uses KSPMonitorLGResidualNormCreate()
1508: . -ksp_monitor_lg_true_residualnorm - sets line graph monitor,
1509: uses KSPMonitorLGResidualNormCreate()
1510: . -ksp_monitor_singular_value - sets KSPMonitorSingularValue()
1511: - -ksp_monitor_cancel - cancels all monitors that have
1512: been hardwired into a code by
1513: calls to KSPMonitorSet(), but
1514: does not cancel those set via
1515: the options database.
1517: Notes:
1518: The default is to do nothing. To print the residual, or preconditioned
1519: residual if KSPSetNormType(ksp,KSP_NORM_PRECONDITIONED) was called, use
1520: KSPMonitorDefault() as the monitoring routine, with a null monitoring
1521: context.
1523: Several different monitoring routines may be set by calling
1524: KSPMonitorSet() multiple times; all will be called in the
1525: order in which they were set.
1527: Fortran notes: Only a single monitor function can be set for each KSP object
1529: Level: beginner
1531: .keywords: KSP, set, monitor
1533: .seealso: KSPMonitorDefault(), KSPMonitorLGResidualNormCreate(), KSPMonitorCancel()
1534: @*/
1535: PetscErrorCode KSPMonitorSet(KSP ksp,PetscErrorCode (*monitor)(KSP,PetscInt,PetscReal,void*),void *mctx,PetscErrorCode (*monitordestroy)(void**))
1536: {
1537: PetscInt i;
1542: if (ksp->numbermonitors >= MAXKSPMONITORS) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_OUTOFRANGE,"Too many KSP monitors set");
1543: for (i=0; i<ksp->numbermonitors;i++) {
1544: if (monitor == ksp->monitor[i] && monitordestroy == ksp->monitordestroy[i] && mctx == ksp->monitorcontext[i]) {
1545: if (monitordestroy) {
1546: (*monitordestroy)(&mctx);
1547: }
1548: return(0);
1549: }
1550: }
1551: ksp->monitor[ksp->numbermonitors] = monitor;
1552: ksp->monitordestroy[ksp->numbermonitors] = monitordestroy;
1553: ksp->monitorcontext[ksp->numbermonitors++] = (void*)mctx;
1554: return(0);
1555: }
1559: /*@
1560: KSPMonitorCancel - Clears all monitors for a KSP object.
1562: Logically Collective on KSP
1564: Input Parameters:
1565: . ksp - iterative context obtained from KSPCreate()
1567: Options Database Key:
1568: . -ksp_monitor_cancel - Cancels all monitors that have
1569: been hardwired into a code by calls to KSPMonitorSet(),
1570: but does not cancel those set via the options database.
1572: Level: intermediate
1574: .keywords: KSP, set, monitor
1576: .seealso: KSPMonitorDefault(), KSPMonitorLGResidualNormCreate(), KSPMonitorSet()
1577: @*/
1578: PetscErrorCode KSPMonitorCancel(KSP ksp)
1579: {
1581: PetscInt i;
1585: for (i=0; i<ksp->numbermonitors; i++) {
1586: if (ksp->monitordestroy[i]) {
1587: (*ksp->monitordestroy[i])(&ksp->monitorcontext[i]);
1588: }
1589: }
1590: ksp->numbermonitors = 0;
1591: return(0);
1592: }
1596: /*@C
1597: KSPGetMonitorContext - Gets the monitoring context, as set by
1598: KSPMonitorSet() for the FIRST monitor only.
1600: Not Collective
1602: Input Parameter:
1603: . ksp - iterative context obtained from KSPCreate()
1605: Output Parameter:
1606: . ctx - monitoring context
1608: Level: intermediate
1610: .keywords: KSP, get, monitor, context
1612: .seealso: KSPMonitorDefault(), KSPMonitorLGResidualNormCreate()
1613: @*/
1614: PetscErrorCode KSPGetMonitorContext(KSP ksp,void **ctx)
1615: {
1618: *ctx = (ksp->monitorcontext[0]);
1619: return(0);
1620: }
1624: /*@
1625: KSPSetResidualHistory - Sets the array used to hold the residual history.
1626: If set, this array will contain the residual norms computed at each
1627: iteration of the solver.
1629: Not Collective
1631: Input Parameters:
1632: + ksp - iterative context obtained from KSPCreate()
1633: . a - array to hold history
1634: . na - size of a
1635: - reset - PETSC_TRUE indicates the history counter is reset to zero
1636: for each new linear solve
1638: Level: advanced
1640: Notes: The array is NOT freed by PETSc so the user needs to keep track of
1641: it and destroy once the KSP object is destroyed.
1643: If 'a' is NULL then space is allocated for the history. If 'na' PETSC_DECIDE or PETSC_DEFAULT then a
1644: default array of length 10000 is allocated.
1646: .keywords: KSP, set, residual, history, norm
1648: .seealso: KSPGetResidualHistory()
1650: @*/
1651: PetscErrorCode KSPSetResidualHistory(KSP ksp,PetscReal a[],PetscInt na,PetscBool reset)
1652: {
1658: PetscFree(ksp->res_hist_alloc);
1659: if (na != PETSC_DECIDE && na != PETSC_DEFAULT && a) {
1660: ksp->res_hist = a;
1661: ksp->res_hist_max = na;
1662: } else {
1663: if (na != PETSC_DECIDE && na != PETSC_DEFAULT) ksp->res_hist_max = na;
1664: else ksp->res_hist_max = 10000; /* like default ksp->max_it */
1665: PetscCalloc1(ksp->res_hist_max,&ksp->res_hist_alloc);
1667: ksp->res_hist = ksp->res_hist_alloc;
1668: }
1669: ksp->res_hist_len = 0;
1670: ksp->res_hist_reset = reset;
1671: return(0);
1672: }
1676: /*@C
1677: KSPGetResidualHistory - Gets the array used to hold the residual history
1678: and the number of residuals it contains.
1680: Not Collective
1682: Input Parameter:
1683: . ksp - iterative context obtained from KSPCreate()
1685: Output Parameters:
1686: + a - pointer to array to hold history (or NULL)
1687: - na - number of used entries in a (or NULL)
1689: Level: advanced
1691: Notes:
1692: Can only be called after a KSPSetResidualHistory() otherwise a and na are set to zero
1694: The Fortran version of this routine has a calling sequence
1695: $ call KSPGetResidualHistory(KSP ksp, integer na, integer ierr)
1696: note that you have passed a Fortran array into KSPSetResidualHistory() and you need
1697: to access the residual values from this Fortran array you provided. Only the na (number of
1698: residual norms currently held) is set.
1700: .keywords: KSP, get, residual, history, norm
1702: .seealso: KSPGetResidualHistory()
1704: @*/
1705: PetscErrorCode KSPGetResidualHistory(KSP ksp,PetscReal *a[],PetscInt *na)
1706: {
1709: if (a) *a = ksp->res_hist;
1710: if (na) *na = ksp->res_hist_len;
1711: return(0);
1712: }
1716: /*@C
1717: KSPSetConvergenceTest - Sets the function to be used to determine
1718: convergence.
1720: Logically Collective on KSP
1722: Input Parameters:
1723: + ksp - iterative context obtained from KSPCreate()
1724: . converge - pointer to int function
1725: . cctx - context for private data for the convergence routine (may be null)
1726: - destroy - a routine for destroying the context (may be null)
1728: Calling sequence of converge:
1729: $ converge (KSP ksp, int it, PetscReal rnorm, KSPConvergedReason *reason,void *mctx)
1731: + ksp - iterative context obtained from KSPCreate()
1732: . it - iteration number
1733: . rnorm - (estimated) 2-norm of (preconditioned) residual
1734: . reason - the reason why it has converged or diverged
1735: - cctx - optional convergence context, as set by KSPSetConvergenceTest()
1738: Notes:
1739: Must be called after the KSP type has been set so put this after
1740: a call to KSPSetType(), or KSPSetFromOptions().
1742: The default convergence test, KSPConvergedDefault(), aborts if the
1743: residual grows to more than 10000 times the initial residual.
1745: The default is a combination of relative and absolute tolerances.
1746: The residual value that is tested may be an approximation; routines
1747: that need exact values should compute them.
1749: In the default PETSc convergence test, the precise values of reason
1750: are macros such as KSP_CONVERGED_RTOL, which are defined in petscksp.h.
1752: Level: advanced
1754: .keywords: KSP, set, convergence, test, context
1756: .seealso: KSPConvergedDefault(), KSPGetConvergenceContext(), KSPSetTolerances()
1757: @*/
1758: PetscErrorCode KSPSetConvergenceTest(KSP ksp,PetscErrorCode (*converge)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*),void *cctx,PetscErrorCode (*destroy)(void*))
1759: {
1764: if (ksp->convergeddestroy) {
1765: (*ksp->convergeddestroy)(ksp->cnvP);
1766: }
1767: ksp->converged = converge;
1768: ksp->convergeddestroy = destroy;
1769: ksp->cnvP = (void*)cctx;
1770: return(0);
1771: }
1775: /*@C
1776: KSPGetConvergenceContext - Gets the convergence context set with
1777: KSPSetConvergenceTest().
1779: Not Collective
1781: Input Parameter:
1782: . ksp - iterative context obtained from KSPCreate()
1784: Output Parameter:
1785: . ctx - monitoring context
1787: Level: advanced
1789: .keywords: KSP, get, convergence, test, context
1791: .seealso: KSPConvergedDefault(), KSPSetConvergenceTest()
1792: @*/
1793: PetscErrorCode KSPGetConvergenceContext(KSP ksp,void **ctx)
1794: {
1797: *ctx = ksp->cnvP;
1798: return(0);
1799: }
1803: /*@C
1804: KSPBuildSolution - Builds the approximate solution in a vector provided.
1805: This routine is NOT commonly needed (see KSPSolve()).
1807: Collective on KSP
1809: Input Parameter:
1810: . ctx - iterative context obtained from KSPCreate()
1812: Output Parameter:
1813: Provide exactly one of
1814: + v - location to stash solution.
1815: - V - the solution is returned in this location. This vector is created
1816: internally. This vector should NOT be destroyed by the user with
1817: VecDestroy().
1819: Notes:
1820: This routine can be used in one of two ways
1821: .vb
1822: KSPBuildSolution(ksp,NULL,&V);
1823: or
1824: KSPBuildSolution(ksp,v,NULL); or KSPBuildSolution(ksp,v,&v);
1825: .ve
1826: In the first case an internal vector is allocated to store the solution
1827: (the user cannot destroy this vector). In the second case the solution
1828: is generated in the vector that the user provides. Note that for certain
1829: methods, such as KSPCG, the second case requires a copy of the solution,
1830: while in the first case the call is essentially free since it simply
1831: returns the vector where the solution already is stored. For some methods
1832: like GMRES this is a reasonably expensive operation and should only be
1833: used in truly needed.
1835: Level: advanced
1837: .keywords: KSP, build, solution
1839: .seealso: KSPGetSolution(), KSPBuildResidual()
1840: @*/
1841: PetscErrorCode KSPBuildSolution(KSP ksp,Vec v,Vec *V)
1842: {
1847: if (!V && !v) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONG,"Must provide either v or V");
1848: if (!V) V = &v;
1849: (*ksp->ops->buildsolution)(ksp,v,V);
1850: return(0);
1851: }
1855: /*@C
1856: KSPBuildResidual - Builds the residual in a vector provided.
1858: Collective on KSP
1860: Input Parameter:
1861: . ksp - iterative context obtained from KSPCreate()
1863: Output Parameters:
1864: + v - optional location to stash residual. If v is not provided,
1865: then a location is generated.
1866: . t - work vector. If not provided then one is generated.
1867: - V - the residual
1869: Notes:
1870: Regardless of whether or not v is provided, the residual is
1871: returned in V.
1873: Level: advanced
1875: .keywords: KSP, build, residual
1877: .seealso: KSPBuildSolution()
1878: @*/
1879: PetscErrorCode KSPBuildResidual(KSP ksp,Vec t,Vec v,Vec *V)
1880: {
1882: PetscBool flag = PETSC_FALSE;
1883: Vec w = v,tt = t;
1887: if (!w) {
1888: VecDuplicate(ksp->vec_rhs,&w);
1889: PetscLogObjectParent((PetscObject)ksp,(PetscObject)w);
1890: }
1891: if (!tt) {
1892: VecDuplicate(ksp->vec_sol,&tt); flag = PETSC_TRUE;
1893: PetscLogObjectParent((PetscObject)ksp,(PetscObject)tt);
1894: }
1895: (*ksp->ops->buildresidual)(ksp,tt,w,V);
1896: if (flag) {VecDestroy(&tt);}
1897: return(0);
1898: }
1902: /*@
1903: KSPSetDiagonalScale - Tells KSP to symmetrically diagonally scale the system
1904: before solving. This actually CHANGES the matrix (and right hand side).
1906: Logically Collective on KSP
1908: Input Parameter:
1909: + ksp - the KSP context
1910: - scale - PETSC_TRUE or PETSC_FALSE
1912: Options Database Key:
1913: + -ksp_diagonal_scale -
1914: - -ksp_diagonal_scale_fix - scale the matrix back AFTER the solve
1917: Notes: Scales the matrix by D^(-1/2) A D^(-1/2) [D^(1/2) x ] = D^(-1/2) b
1918: where D_{ii} is 1/abs(A_{ii}) unless A_{ii} is zero and then it is 1.
1920: BE CAREFUL with this routine: it actually scales the matrix and right
1921: hand side that define the system. After the system is solved the matrix
1922: and right hand side remain scaled unless you use KSPSetDiagonalScaleFix()
1924: This should NOT be used within the SNES solves if you are using a line
1925: search.
1927: If you use this with the PCType Eisenstat preconditioner than you can
1928: use the PCEisenstatNoDiagonalScaling() option, or -pc_eisenstat_no_diagonal_scaling
1929: to save some unneeded, redundant flops.
1931: Level: intermediate
1933: .keywords: KSP, set, options, prefix, database
1935: .seealso: KSPGetDiagonalScale(), KSPSetDiagonalScaleFix()
1936: @*/
1937: PetscErrorCode KSPSetDiagonalScale(KSP ksp,PetscBool scale)
1938: {
1942: ksp->dscale = scale;
1943: return(0);
1944: }
1948: /*@
1949: KSPGetDiagonalScale - Checks if KSP solver scales the matrix and
1950: right hand side
1952: Not Collective
1954: Input Parameter:
1955: . ksp - the KSP context
1957: Output Parameter:
1958: . scale - PETSC_TRUE or PETSC_FALSE
1960: Notes:
1961: BE CAREFUL with this routine: it actually scales the matrix and right
1962: hand side that define the system. After the system is solved the matrix
1963: and right hand side remain scaled unless you use KSPSetDiagonalScaleFix()
1965: Level: intermediate
1967: .keywords: KSP, set, options, prefix, database
1969: .seealso: KSPSetDiagonalScale(), KSPSetDiagonalScaleFix()
1970: @*/
1971: PetscErrorCode KSPGetDiagonalScale(KSP ksp,PetscBool *scale)
1972: {
1976: *scale = ksp->dscale;
1977: return(0);
1978: }
1982: /*@
1983: KSPSetDiagonalScaleFix - Tells KSP to diagonally scale the system
1984: back after solving.
1986: Logically Collective on KSP
1988: Input Parameter:
1989: + ksp - the KSP context
1990: - fix - PETSC_TRUE to scale back after the system solve, PETSC_FALSE to not
1991: rescale (default)
1993: Notes:
1994: Must be called after KSPSetDiagonalScale()
1996: Using this will slow things down, because it rescales the matrix before and
1997: after each linear solve. This is intended mainly for testing to allow one
1998: to easily get back the original system to make sure the solution computed is
1999: accurate enough.
2001: Level: intermediate
2003: .keywords: KSP, set, options, prefix, database
2005: .seealso: KSPGetDiagonalScale(), KSPSetDiagonalScale(), KSPGetDiagonalScaleFix()
2006: @*/
2007: PetscErrorCode KSPSetDiagonalScaleFix(KSP ksp,PetscBool fix)
2008: {
2012: ksp->dscalefix = fix;
2013: return(0);
2014: }
2018: /*@
2019: KSPGetDiagonalScaleFix - Determines if KSP diagonally scales the system
2020: back after solving.
2022: Not Collective
2024: Input Parameter:
2025: . ksp - the KSP context
2027: Output Parameter:
2028: . fix - PETSC_TRUE to scale back after the system solve, PETSC_FALSE to not
2029: rescale (default)
2031: Notes:
2032: Must be called after KSPSetDiagonalScale()
2034: If PETSC_TRUE will slow things down, because it rescales the matrix before and
2035: after each linear solve. This is intended mainly for testing to allow one
2036: to easily get back the original system to make sure the solution computed is
2037: accurate enough.
2039: Level: intermediate
2041: .keywords: KSP, set, options, prefix, database
2043: .seealso: KSPGetDiagonalScale(), KSPSetDiagonalScale(), KSPSetDiagonalScaleFix()
2044: @*/
2045: PetscErrorCode KSPGetDiagonalScaleFix(KSP ksp,PetscBool *fix)
2046: {
2050: *fix = ksp->dscalefix;
2051: return(0);
2052: }
2056: /*@C
2057: KSPSetComputeOperators - set routine to compute the linear operators
2059: Logically Collective
2061: Input Arguments:
2062: + ksp - the KSP context
2063: . func - function to compute the operators
2064: - ctx - optional context
2066: Calling sequence of func:
2067: $ func(KSP ksp,Mat A,Mat B,void *ctx)
2069: + ksp - the KSP context
2070: . A - the linear operator
2071: . B - preconditioning matrix
2072: - ctx - optional user-provided context
2074: Notes: The user provided func() will be called automatically at the very next call to KSPSolve(). It will not be called at future KSPSolve() calls
2075: unless either KSPSetComputeOperators() or KSPSetOperators() is called before that KSPSolve() is called.
2077: To reuse the same preconditioner for the next KSPSolve() and not compute a new one based on the most recently computed matrix call KSPSetReusePreconditioner()
2079: Level: beginner
2081: .seealso: KSPSetOperators(), KSPSetComputeRHS(), DMKSPSetComputeOperators(), KSPSetComputeInitialGuess()
2082: @*/
2083: PetscErrorCode KSPSetComputeOperators(KSP ksp,PetscErrorCode (*func)(KSP,Mat,Mat,void*),void *ctx)
2084: {
2086: DM dm;
2090: KSPGetDM(ksp,&dm);
2091: DMKSPSetComputeOperators(dm,func,ctx);
2092: if (ksp->setupstage == KSP_SETUP_NEWRHS) ksp->setupstage = KSP_SETUP_NEWMATRIX;
2093: return(0);
2094: }
2098: /*@C
2099: KSPSetComputeRHS - set routine to compute the right hand side of the linear system
2101: Logically Collective
2103: Input Arguments:
2104: + ksp - the KSP context
2105: . func - function to compute the right hand side
2106: - ctx - optional context
2108: Calling sequence of func:
2109: $ func(KSP ksp,Vec b,void *ctx)
2111: + ksp - the KSP context
2112: . b - right hand side of linear system
2113: - ctx - optional user-provided context
2115: Notes: The routine you provide will be called EACH you call KSPSolve() to prepare the new right hand side for that solve
2117: Level: beginner
2119: .seealso: KSPSolve(), DMKSPSetComputeRHS(), KSPSetComputeOperators()
2120: @*/
2121: PetscErrorCode KSPSetComputeRHS(KSP ksp,PetscErrorCode (*func)(KSP,Vec,void*),void *ctx)
2122: {
2124: DM dm;
2128: KSPGetDM(ksp,&dm);
2129: DMKSPSetComputeRHS(dm,func,ctx);
2130: return(0);
2131: }
2135: /*@C
2136: KSPSetComputeInitialGuess - set routine to compute the initial guess of the linear system
2138: Logically Collective
2140: Input Arguments:
2141: + ksp - the KSP context
2142: . func - function to compute the initial guess
2143: - ctx - optional context
2145: Calling sequence of func:
2146: $ func(KSP ksp,Vec x,void *ctx)
2148: + ksp - the KSP context
2149: . x - solution vector
2150: - ctx - optional user-provided context
2152: Level: beginner
2154: .seealso: KSPSolve(), KSPSetComputeRHS(), KSPSetComputeOperators(), DMKSPSetComputeInitialGuess()
2155: @*/
2156: PetscErrorCode KSPSetComputeInitialGuess(KSP ksp,PetscErrorCode (*func)(KSP,Vec,void*),void *ctx)
2157: {
2159: DM dm;
2163: KSPGetDM(ksp,&dm);
2164: DMKSPSetComputeInitialGuess(dm,func,ctx);
2165: return(0);
2166: }