Actual source code: cheby.c
petsc-3.5.4 2015-05-23
2: #include <petsc-private/kspimpl.h> /*I "petscksp.h" I*/
3: #include <../src/ksp/ksp/impls/cheby/chebyshevimpl.h>
7: PetscErrorCode KSPReset_Chebyshev(KSP ksp)
8: {
9: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
13: KSPReset(cheb->kspest);
14: return(0);
15: }
19: PetscErrorCode KSPSetUp_Chebyshev(KSP ksp)
20: {
21: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
25: KSPSetWorkVecs(ksp,3);
26: if (cheb->emin == 0. || cheb->emax == 0.) { /* We need to estimate eigenvalues */
27: KSPChebyshevSetEstimateEigenvalues(ksp,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE);
28: /* Enable runtime options for cheb->kspest:
29: KSPChebyshevSetEstimateEigenvalues() creates cheb->kspest, but does not call KSPSetFromOptions(cheb->kspest)! */
30: KSPSetFromOptions(cheb->kspest);
31: }
32: return(0);
33: }
37: static PetscErrorCode KSPChebyshevSetEigenvalues_Chebyshev(KSP ksp,PetscReal emax,PetscReal emin)
38: {
39: KSP_Chebyshev *chebyshevP = (KSP_Chebyshev*)ksp->data;
43: if (emax <= emin) SETERRQ2(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"Maximum eigenvalue must be larger than minimum: max %g min %g",(double)emax,(double)emin);
44: if (emax*emin <= 0.0) SETERRQ2(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"Both eigenvalues must be of the same sign: max %g min %g",(double)emax,(double)emin);
45: chebyshevP->emax = emax;
46: chebyshevP->emin = emin;
48: KSPChebyshevSetEstimateEigenvalues(ksp,0.,0.,0.,0.); /* Destroy any estimation setup */
49: return(0);
50: }
54: static PetscErrorCode KSPChebyshevSetEstimateEigenvalues_Chebyshev(KSP ksp,PetscReal a,PetscReal b,PetscReal c,PetscReal d)
55: {
56: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
60: if (a != 0.0 || b != 0.0 || c != 0.0 || d != 0.0) {
61: if (!cheb->kspest) { /* should this block of code be moved to KSPSetUp_Chebyshev()? */
62: PetscBool nonzero;
64: KSPCreate(PetscObjectComm((PetscObject)ksp),&cheb->kspest);
65: PetscObjectIncrementTabLevel((PetscObject)cheb->kspest,(PetscObject)ksp,1);
66: KSPSetOptionsPrefix(cheb->kspest,((PetscObject)ksp)->prefix);
67: KSPAppendOptionsPrefix(cheb->kspest,"est_");
69: KSPGetPC(cheb->kspest,&cheb->pcnone);
70: PetscObjectReference((PetscObject)cheb->pcnone);
71: PCSetType(cheb->pcnone,PCNONE);
72: KSPSetPC(cheb->kspest,ksp->pc);
73:
74: KSPGetInitialGuessNonzero(ksp,&nonzero);
75: KSPSetInitialGuessNonzero(cheb->kspest,nonzero);
76: KSPSetComputeEigenvalues(cheb->kspest,PETSC_TRUE);
78: /* Estimate with a fixed number of iterations */
79: KSPSetConvergenceTest(cheb->kspest,KSPConvergedSkip,0,0);
80: KSPSetNormType(cheb->kspest,KSP_NORM_NONE);
81: KSPSetTolerances(cheb->kspest,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,cheb->eststeps);
82: }
83: if (a >= 0) cheb->tform[0] = a;
84: if (b >= 0) cheb->tform[1] = b;
85: if (c >= 0) cheb->tform[2] = c;
86: if (d >= 0) cheb->tform[3] = d;
87: cheb->estimate_current = PETSC_FALSE;
88: } else {
89: KSPDestroy(&cheb->kspest);
90: PCDestroy(&cheb->pcnone);
91: }
92: return(0);
93: }
97: static PetscErrorCode KSPChebyshevEstEigSetRandom_Chebyshev(KSP ksp,PetscRandom random)
98: {
99: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
103: if (random) {PetscObjectReference((PetscObject)random);}
104: PetscRandomDestroy(&cheb->random);
106: cheb->random = random;
107: return(0);
108: }
112: static PetscErrorCode KSPChebyshevSetNewMatrix_Chebyshev(KSP ksp)
113: {
114: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
117: cheb->estimate_current = PETSC_FALSE;
118: return(0);
119: }
123: /*@
124: KSPChebyshevSetEigenvalues - Sets estimates for the extreme eigenvalues
125: of the preconditioned problem.
127: Logically Collective on KSP
129: Input Parameters:
130: + ksp - the Krylov space context
131: - emax, emin - the eigenvalue estimates
133: Options Database:
134: . -ksp_chebyshev_eigenvalues emin,emax
136: Note: If you run with the Krylov method of KSPCG with the option -ksp_monitor_singular_value it will
137: for that given matrix and preconditioner an estimate of the extreme eigenvalues.
139: Level: intermediate
141: .keywords: KSP, Chebyshev, set, eigenvalues
142: @*/
143: PetscErrorCode KSPChebyshevSetEigenvalues(KSP ksp,PetscReal emax,PetscReal emin)
144: {
151: PetscTryMethod(ksp,"KSPChebyshevSetEigenvalues_C",(KSP,PetscReal,PetscReal),(ksp,emax,emin));
152: return(0);
153: }
157: /*@
158: KSPChebyshevSetEstimateEigenvalues - Automatically estimate the eigenvalues to use for Chebyshev
160: Logically Collective on KSP
162: Input Parameters:
163: + ksp - the Krylov space context
164: . a - multiple of min eigenvalue estimate to use for min Chebyshev bound (or PETSC_DECIDE)
165: . b - multiple of max eigenvalue estimate to use for min Chebyshev bound (or PETSC_DECIDE)
166: . c - multiple of min eigenvalue estimate to use for max Chebyshev bound (or PETSC_DECIDE)
167: - d - multiple of max eigenvalue estimate to use for max Chebyshev bound (or PETSC_DECIDE)
169: Options Database:
170: . -ksp_chebyshev_estimate_eigenvalues a,b,c,d
172: Notes:
173: The Chebyshev bounds are estimated using
174: .vb
175: minbound = a*minest + b*maxest
176: maxbound = c*minest + d*maxest
177: .ve
178: The default configuration targets the upper part of the spectrum for use as a multigrid smoother, so only the maximum eigenvalue estimate is used.
179: The minimum eigenvalue estimate obtained by Krylov iteration is typically not accurate until the method has converged.
181: If 0.0 is passed for all transform arguments (a,b,c,d), eigenvalue estimation is disabled.
183: The default transform is (0,0.1; 0,1.1) which targets the "upper" part of the spectrum, as desirable for use with multigrid.
185: Level: intermediate
187: .keywords: KSP, Chebyshev, set, eigenvalues, PCMG
188: @*/
189: PetscErrorCode KSPChebyshevSetEstimateEigenvalues(KSP ksp,PetscReal a,PetscReal b,PetscReal c,PetscReal d)
190: {
199: PetscTryMethod(ksp,"KSPChebyshevSetEstimateEigenvalues_C",(KSP,PetscReal,PetscReal,PetscReal,PetscReal),(ksp,a,b,c,d));
200: return(0);
201: }
205: /*@
206: KSPChebyshevEstEigSetRandom - set random context for estimating eigenvalues
208: Logically Collective
210: Input Arguments:
211: + ksp - linear solver context
212: - random - random number context or NULL to disable randomized RHS
214: Level: intermediate
216: .seealso: KSPChebyshevSetEstimateEigenvalues(), PetscRandomCreate()
217: @*/
218: PetscErrorCode KSPChebyshevEstEigSetRandom(KSP ksp,PetscRandom random)
219: {
225: PetscTryMethod(ksp,"KSPChebyshevEstEigSetRandom_C",(KSP,PetscRandom),(ksp,random));
226: return(0);
227: }
231: /*@
232: KSPChebyshevSetNewMatrix - Indicates that the matrix has changed, causes eigenvalue estimates to be recomputed if appropriate.
234: Logically Collective on KSP
236: Input Parameter:
237: . ksp - the Krylov space context
239: Level: developer
241: .keywords: KSP, Chebyshev, set, eigenvalues
243: .seealso: KSPChebyshevSetEstimateEigenvalues()
244: @*/
245: PetscErrorCode KSPChebyshevSetNewMatrix(KSP ksp)
246: {
251: PetscTryMethod(ksp,"KSPChebyshevSetNewMatrix_C",(KSP),(ksp));
252: return(0);
253: }
257: PetscErrorCode KSPSetFromOptions_Chebyshev(KSP ksp)
258: {
259: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
261: PetscInt two = 2,four = 4;
262: PetscReal tform[4] = {PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE};
263: PetscBool flg;
266: PetscOptionsHead("KSP Chebyshev Options");
267: PetscOptionsInt("-ksp_chebyshev_eststeps","Number of est steps in Chebyshev","",cheb->eststeps,&cheb->eststeps,NULL);
268: PetscOptionsRealArray("-ksp_chebyshev_eigenvalues","extreme eigenvalues","KSPChebyshevSetEigenvalues",&cheb->emin,&two,0);
269: PetscOptionsRealArray("-ksp_chebyshev_estimate_eigenvalues","estimate eigenvalues using a Krylov method, then use this transform for Chebyshev eigenvalue bounds","KSPChebyshevSetEstimateEigenvalues",tform,&four,&flg);
270: if (flg) {
271: switch (four) {
272: case 0:
273: KSPChebyshevSetEstimateEigenvalues(ksp,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE);
274: break;
275: case 2: /* Base everything on the max eigenvalues */
276: KSPChebyshevSetEstimateEigenvalues(ksp,PETSC_DECIDE,tform[0],PETSC_DECIDE,tform[1]);
277: break;
278: case 4: /* Use the full 2x2 linear transformation */
279: KSPChebyshevSetEstimateEigenvalues(ksp,tform[0],tform[1],tform[2],tform[3]);
280: break;
281: default: SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"Must specify either 0, 2, or 4 parameters for eigenvalue estimation");
282: }
283: }
285: if (cheb->kspest) {
286: PetscBool estrand = PETSC_FALSE;
287: PetscOptionsBool("-ksp_chebyshev_estimate_eigenvalues_random","Use Random right hand side for eigenvalue estimation","KSPChebyshevEstEigSetRandom",estrand,&estrand,NULL);
288: if (estrand) {
289: PetscRandom random;
290: PetscRandomCreate(PetscObjectComm((PetscObject)ksp),&random);
291: PetscObjectSetOptionsPrefix((PetscObject)random,((PetscObject)ksp)->prefix);
292: PetscObjectAppendOptionsPrefix((PetscObject)random,"ksp_chebyshev_estimate_eigenvalues_");
293: PetscRandomSetFromOptions(random);
294: KSPChebyshevEstEigSetRandom(ksp,random);
295: PetscRandomDestroy(&random);
296: }
297: }
299: if (cheb->kspest) {
300: /* Mask the PC so that PCSetFromOptions does not do anything */
301: KSPSetPC(cheb->kspest,cheb->pcnone);
302: KSPSetOptionsPrefix(cheb->kspest,((PetscObject)ksp)->prefix);
303: KSPAppendOptionsPrefix(cheb->kspest,"est_");
304: if (!((PetscObject)cheb->kspest)->type_name) {
305: KSPSetType(cheb->kspest,KSPGMRES);
306: }
307: KSPSetFromOptions(cheb->kspest);
308: KSPSetPC(cheb->kspest,ksp->pc);
309: }
310: PetscOptionsTail();
311: return(0);
312: }
316: /*
317: * Must be passed a KSP solver that has "converged", with KSPSetComputeEigenvalues() called before the solve
318: */
319: static PetscErrorCode KSPChebyshevComputeExtremeEigenvalues_Private(KSP kspest,PetscReal *emin,PetscReal *emax)
320: {
322: PetscInt n,neig;
323: PetscReal *re,*im,min,max;
326: KSPGetIterationNumber(kspest,&n);
327: PetscMalloc2(n,&re,n,&im);
328: KSPComputeEigenvalues(kspest,n,re,im,&neig);
329: min = PETSC_MAX_REAL;
330: max = PETSC_MIN_REAL;
331: for (n=0; n<neig; n++) {
332: min = PetscMin(min,re[n]);
333: max = PetscMax(max,re[n]);
334: }
335: PetscFree2(re,im);
336: *emax = max;
337: *emin = min;
338: return(0);
339: }
343: PetscErrorCode KSPSolve_Chebyshev(KSP ksp)
344: {
345: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
347: PetscInt k,kp1,km1,maxit,ktmp,i;
348: PetscScalar alpha,omegaprod,mu,omega,Gamma,c[3],scale;
349: PetscReal rnorm = 0.0;
350: Vec sol_orig,b,p[3],r;
351: Mat Amat,Pmat;
352: PetscBool diagonalscale;
355: PCGetDiagonalScale(ksp->pc,&diagonalscale);
356: if (diagonalscale) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Krylov method %s does not support diagonal scaling",((PetscObject)ksp)->type_name);
358: if (cheb->kspest && !cheb->estimate_current) {
359: PetscReal max=0.0,min=0.0;
360: Vec X,B;
361: X = ksp->work[0];
362: if (cheb->random) {
363: B = ksp->work[1];
364: VecSetRandom(B,cheb->random);
365: } else {
366: B = ksp->vec_rhs;
367: }
368: KSPSolve(cheb->kspest,B,X);
369:
370: if (ksp->guess_zero) {
371: VecZeroEntries(X);
372: }
373: KSPChebyshevComputeExtremeEigenvalues_Private(cheb->kspest,&min,&max);
375: cheb->emin = cheb->tform[0]*min + cheb->tform[1]*max;
376: cheb->emax = cheb->tform[2]*min + cheb->tform[3]*max;
378: cheb->estimate_current = PETSC_TRUE;
379: }
381: ksp->its = 0;
382: PCGetOperators(ksp->pc,&Amat,&Pmat);
383: maxit = ksp->max_it;
385: /* These three point to the three active solutions, we
386: rotate these three at each solution update */
387: km1 = 0; k = 1; kp1 = 2;
388: sol_orig = ksp->vec_sol; /* ksp->vec_sol will be asigned to rotating vector p[k], thus save its address */
389: b = ksp->vec_rhs;
390: p[km1] = sol_orig;
391: p[k] = ksp->work[0];
392: p[kp1] = ksp->work[1];
393: r = ksp->work[2];
395: /* use scale*B as our preconditioner */
396: scale = 2.0/(cheb->emax + cheb->emin);
398: /* -alpha <= scale*lambda(B^{-1}A) <= alpha */
399: alpha = 1.0 - scale*(cheb->emin);
400: Gamma = 1.0;
401: mu = 1.0/alpha;
402: omegaprod = 2.0/alpha;
404: c[km1] = 1.0;
405: c[k] = mu;
407: if (!ksp->guess_zero) {
408: KSP_MatMult(ksp,Amat,p[km1],r); /* r = b - A*p[km1] */
409: VecAYPX(r,-1.0,b);
410: } else {
411: VecCopy(b,r);
412: }
414: KSP_PCApply(ksp,r,p[k]); /* p[k] = scale B^{-1}r + p[km1] */
415: VecAYPX(p[k],scale,p[km1]);
417: for (i=0; i<maxit; i++) {
418: PetscObjectSAWsTakeAccess((PetscObject)ksp);
420: ksp->its++;
421: PetscObjectSAWsGrantAccess((PetscObject)ksp);
422: c[kp1] = 2.0*mu*c[k] - c[km1];
423: omega = omegaprod*c[k]/c[kp1];
425: KSP_MatMult(ksp,Amat,p[k],r); /* r = b - Ap[k] */
426: VecAYPX(r,-1.0,b);
427: KSP_PCApply(ksp,r,p[kp1]); /* p[kp1] = B^{-1}r */
428: ksp->vec_sol = p[k];
430: /* calculate residual norm if requested */
431: if (ksp->normtype != KSP_NORM_NONE || ksp->numbermonitors) {
432: if (ksp->normtype == KSP_NORM_UNPRECONDITIONED) {
433: VecNorm(r,NORM_2,&rnorm);
434: } else {
435: VecNorm(p[kp1],NORM_2,&rnorm);
436: }
437: PetscObjectSAWsTakeAccess((PetscObject)ksp);
438: ksp->rnorm = rnorm;
439: PetscObjectSAWsGrantAccess((PetscObject)ksp);
440: KSPLogResidualHistory(ksp,rnorm);
441: KSPMonitor(ksp,i,rnorm);
442: (*ksp->converged)(ksp,i,rnorm,&ksp->reason,ksp->cnvP);
443: if (ksp->reason) break;
444: }
446: /* y^{k+1} = omega(y^{k} - y^{k-1} + Gamma*r^{k}) + y^{k-1} */
447: VecAXPBYPCZ(p[kp1],1.0-omega,omega,omega*Gamma*scale,p[km1],p[k]);
449: ktmp = km1;
450: km1 = k;
451: k = kp1;
452: kp1 = ktmp;
453: }
454: if (!ksp->reason) {
455: if (ksp->normtype != KSP_NORM_NONE) {
456: KSP_MatMult(ksp,Amat,p[k],r); /* r = b - Ap[k] */
457: VecAYPX(r,-1.0,b);
458: if (ksp->normtype == KSP_NORM_UNPRECONDITIONED) {
459: VecNorm(r,NORM_2,&rnorm);
460: } else {
461: KSP_PCApply(ksp,r,p[kp1]); /* p[kp1] = B^{-1}r */
462: VecNorm(p[kp1],NORM_2,&rnorm);
463: }
464: PetscObjectSAWsTakeAccess((PetscObject)ksp);
465: ksp->rnorm = rnorm;
466: PetscObjectSAWsGrantAccess((PetscObject)ksp);
467: ksp->vec_sol = p[k];
468: KSPLogResidualHistory(ksp,rnorm);
469: KSPMonitor(ksp,i,rnorm);
470: }
471: if (ksp->its >= ksp->max_it) {
472: if (ksp->normtype != KSP_NORM_NONE) {
473: (*ksp->converged)(ksp,i,rnorm,&ksp->reason,ksp->cnvP);
474: if (!ksp->reason) ksp->reason = KSP_DIVERGED_ITS;
475: } else ksp->reason = KSP_CONVERGED_ITS;
476: }
477: }
479: /* make sure solution is in vector x */
480: ksp->vec_sol = sol_orig;
481: if (k) {
482: VecCopy(p[k],sol_orig);
483: }
484: return(0);
485: }
489: PetscErrorCode KSPView_Chebyshev(KSP ksp,PetscViewer viewer)
490: {
491: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
493: PetscBool iascii;
496: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
497: if (iascii) {
498: PetscViewerASCIIPrintf(viewer," Chebyshev: eigenvalue estimates: min = %g, max = %g\n",(double)cheb->emin,(double)cheb->emax);
499: if (cheb->kspest) {
500: PetscViewerASCIIPrintf(viewer," Chebyshev: estimated using: [%g %g; %g %g]\n",(double)cheb->tform[0],(double)cheb->tform[1],(double)cheb->tform[2],(double)cheb->tform[3]);
502: if (cheb->random) {
503: PetscViewerASCIIPrintf(viewer," Chebyshev: estimating eigenvalues using random right hand side\n");
504: PetscViewerASCIIPushTab(viewer);
505: PetscRandomView(cheb->random,viewer);
506: PetscViewerASCIIPopTab(viewer);
507: }
508: PetscViewerASCIIPushTab(viewer);
509: KSPView(cheb->kspest,viewer);
510: PetscViewerASCIIPopTab(viewer);
511: }
512: }
513: return(0);
514: }
518: PetscErrorCode KSPDestroy_Chebyshev(KSP ksp)
519: {
520: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
524: KSPDestroy(&cheb->kspest);
525: PCDestroy(&cheb->pcnone);
526: PetscRandomDestroy(&cheb->random);
527: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetEigenvalues_C",NULL);
528: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetEstimateEigenvalues_C",NULL);
529: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevEstEigSetRandom_C",NULL);
530: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetNewMatrix_C",NULL);
531: KSPDestroyDefault(ksp);
532: return(0);
533: }
535: /*MC
536: KSPCHEBYSHEV - The preconditioned Chebyshev iterative method
538: Options Database Keys:
539: + -ksp_chebyshev_eigenvalues <emin,emax> - set approximations to the smallest and largest eigenvalues
540: of the preconditioned operator. If these are accurate you will get much faster convergence.
541: - -ksp_chebyshev_estimate_eigenvalues <a,b,c,d> - estimate eigenvalues using a Krylov method, then use this
542: transform for Chebyshev eigenvalue bounds (KSPChebyshevSetEstimateEigenvalues)
545: Level: beginner
547: Notes: The Chebyshev method requires both the matrix and preconditioner to
548: be symmetric positive (semi) definite.
549: Only support for left preconditioning.
551: Chebyshev is configured as a smoother by default, targetting the "upper" part of the spectrum.
552: The user should call KSPChebyshevSetEigenvalues() if they have eigenvalue estimates.
554: .seealso: KSPCreate(), KSPSetType(), KSPType (for list of available types), KSP,
555: KSPChebyshevSetEigenvalues(), KSPChebyshevSetEstimateEigenvalues(), KSPRICHARDSON, KSPCG, PCMG
557: M*/
561: PETSC_EXTERN PetscErrorCode KSPCreate_Chebyshev(KSP ksp)
562: {
564: KSP_Chebyshev *chebyshevP;
567: PetscNewLog(ksp,&chebyshevP);
569: ksp->data = (void*)chebyshevP;
570: KSPSetSupportedNorm(ksp,KSP_NORM_PRECONDITIONED,PC_LEFT,3);
571: KSPSetSupportedNorm(ksp,KSP_NORM_UNPRECONDITIONED,PC_LEFT,2);
573: chebyshevP->emin = 0.;
574: chebyshevP->emax = 0.;
576: chebyshevP->tform[0] = 0.0;
577: chebyshevP->tform[1] = 0.1;
578: chebyshevP->tform[2] = 0;
579: chebyshevP->tform[3] = 1.1;
580: chebyshevP->eststeps = 10;
581:
582: ksp->ops->setup = KSPSetUp_Chebyshev;
583: ksp->ops->solve = KSPSolve_Chebyshev;
584: ksp->ops->destroy = KSPDestroy_Chebyshev;
585: ksp->ops->buildsolution = KSPBuildSolutionDefault;
586: ksp->ops->buildresidual = KSPBuildResidualDefault;
587: ksp->ops->setfromoptions = KSPSetFromOptions_Chebyshev;
588: ksp->ops->view = KSPView_Chebyshev;
589: ksp->ops->reset = KSPReset_Chebyshev;
591: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetEigenvalues_C",KSPChebyshevSetEigenvalues_Chebyshev);
592: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetEstimateEigenvalues_C",KSPChebyshevSetEstimateEigenvalues_Chebyshev);
593: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevEstEigSetRandom_C",KSPChebyshevEstEigSetRandom_Chebyshev);
594: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetNewMatrix_C",KSPChebyshevSetNewMatrix_Chebyshev);
595: return(0);
596: }