Actual source code: cheby.c
petsc-3.6.4 2016-04-12
2: #include <petsc/private/kspimpl.h> /*I "petscksp.h" I*/
3: #include <../src/ksp/ksp/impls/cheby/chebyshevimpl.h>
7: static PetscErrorCode KSPReset_Chebyshev(KSP ksp)
8: {
9: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
13: KSPReset(cheb->kspest);
14: return(0);
15: }
19: static 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.) && !cheb->kspest) { /* We need to estimate eigenvalues */
27: KSPChebyshevEstEigSet(ksp,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE);
28: }
29: return(0);
30: }
34: static PetscErrorCode KSPChebyshevSetEigenvalues_Chebyshev(KSP ksp,PetscReal emax,PetscReal emin)
35: {
36: KSP_Chebyshev *chebyshevP = (KSP_Chebyshev*)ksp->data;
40: 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);
41: 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);
42: chebyshevP->emax = emax;
43: chebyshevP->emin = emin;
45: KSPChebyshevEstEigSet(ksp,0.,0.,0.,0.); /* Destroy any estimation setup */
46: return(0);
47: }
51: static PetscErrorCode KSPChebyshevEstEigSet_Chebyshev(KSP ksp,PetscReal a,PetscReal b,PetscReal c,PetscReal d)
52: {
53: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
57: if (a != 0.0 || b != 0.0 || c != 0.0 || d != 0.0) {
58: if (!cheb->kspest) { /* should this block of code be moved to KSPSetUp_Chebyshev()? */
59: PetscBool nonzero;
61: KSPCreate(PetscObjectComm((PetscObject)ksp),&cheb->kspest);
62: PetscObjectIncrementTabLevel((PetscObject)cheb->kspest,(PetscObject)ksp,1);
63: KSPSetOptionsPrefix(cheb->kspest,((PetscObject)ksp)->prefix);
64: KSPAppendOptionsPrefix(cheb->kspest,"esteig_");
65: KSPSetSkipPCSetFromOptions(cheb->kspest,PETSC_TRUE);
67: KSPSetPC(cheb->kspest,ksp->pc);
69: KSPGetInitialGuessNonzero(ksp,&nonzero);
70: KSPSetInitialGuessNonzero(cheb->kspest,nonzero);
71: KSPSetComputeEigenvalues(cheb->kspest,PETSC_TRUE);
73: /* Estimate with a fixed number of iterations */
74: KSPSetConvergenceTest(cheb->kspest,KSPConvergedSkip,0,0);
75: KSPSetNormType(cheb->kspest,KSP_NORM_NONE);
76: KSPSetTolerances(cheb->kspest,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,cheb->eststeps);
77: }
78: if (a >= 0) cheb->tform[0] = a;
79: if (b >= 0) cheb->tform[1] = b;
80: if (c >= 0) cheb->tform[2] = c;
81: if (d >= 0) cheb->tform[3] = d;
82: cheb->amatid = 0;
83: cheb->pmatid = 0;
84: cheb->amatstate = -1;
85: cheb->pmatstate = -1;
86: } else {
87: KSPDestroy(&cheb->kspest);
88: }
89: return(0);
90: }
94: static PetscErrorCode KSPChebyshevEstEigSetRandom_Chebyshev(KSP ksp,PetscRandom random)
95: {
96: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
100: if (random) {PetscObjectReference((PetscObject)random);}
101: PetscRandomDestroy(&cheb->random);
103: cheb->random = random;
104: return(0);
105: }
109: /*@
110: KSPChebyshevSetEigenvalues - Sets estimates for the extreme eigenvalues
111: of the preconditioned problem.
113: Logically Collective on KSP
115: Input Parameters:
116: + ksp - the Krylov space context
117: - emax, emin - the eigenvalue estimates
119: Options Database:
120: . -ksp_chebyshev_eigenvalues emin,emax
122: Note: If you run with the Krylov method of KSPCG with the option -ksp_monitor_singular_value it will
123: for that given matrix and preconditioner an estimate of the extreme eigenvalues.
125: Level: intermediate
127: .keywords: KSP, Chebyshev, set, eigenvalues
128: @*/
129: PetscErrorCode KSPChebyshevSetEigenvalues(KSP ksp,PetscReal emax,PetscReal emin)
130: {
137: PetscTryMethod(ksp,"KSPChebyshevSetEigenvalues_C",(KSP,PetscReal,PetscReal),(ksp,emax,emin));
138: return(0);
139: }
143: /*@
144: KSPChebyshevEstEigSet - Automatically estimate the eigenvalues to use for Chebyshev
146: Logically Collective on KSP
148: Input Parameters:
149: + ksp - the Krylov space context
150: . a - multiple of min eigenvalue estimate to use for min Chebyshev bound (or PETSC_DECIDE)
151: . b - multiple of max eigenvalue estimate to use for min Chebyshev bound (or PETSC_DECIDE)
152: . c - multiple of min eigenvalue estimate to use for max Chebyshev bound (or PETSC_DECIDE)
153: - d - multiple of max eigenvalue estimate to use for max Chebyshev bound (or PETSC_DECIDE)
155: Options Database:
156: . -ksp_chebyshev_esteig a,b,c,d
158: Notes:
159: The Chebyshev bounds are estimated using
160: .vb
161: minbound = a*minest + b*maxest
162: maxbound = c*minest + d*maxest
163: .ve
164: The default configuration targets the upper part of the spectrum for use as a multigrid smoother, so only the maximum eigenvalue estimate is used.
165: The minimum eigenvalue estimate obtained by Krylov iteration is typically not accurate until the method has converged.
167: If 0.0 is passed for all transform arguments (a,b,c,d), eigenvalue estimation is disabled.
169: The default transform is (0,0.1; 0,1.1) which targets the "upper" part of the spectrum, as desirable for use with multigrid.
171: Level: intermediate
173: .keywords: KSP, Chebyshev, set, eigenvalues, PCMG
174: @*/
175: PetscErrorCode KSPChebyshevEstEigSet(KSP ksp,PetscReal a,PetscReal b,PetscReal c,PetscReal d)
176: {
185: PetscTryMethod(ksp,"KSPChebyshevEstEigSet_C",(KSP,PetscReal,PetscReal,PetscReal,PetscReal),(ksp,a,b,c,d));
186: return(0);
187: }
191: /*@
192: KSPChebyshevEstEigSetRandom - set random context for estimating eigenvalues
194: Logically Collective
196: Input Arguments:
197: + ksp - linear solver context
198: - random - random number context or NULL to disable randomized RHS
200: Options Database:
201: . -ksp_chebyshev_esteig_random
203: Level: intermediate
205: .seealso: KSPChebyshevEstEigSet(), PetscRandomCreate()
206: @*/
207: PetscErrorCode KSPChebyshevEstEigSetRandom(KSP ksp,PetscRandom random)
208: {
214: PetscTryMethod(ksp,"KSPChebyshevEstEigSetRandom_C",(KSP,PetscRandom),(ksp,random));
215: return(0);
216: }
220: /*@
221: KSPChebyshevEstEigGetKSP - Get the Krylov method context used to estimate eigenvalues for the Chebyshev method. If
222: a Krylov method is not being used for this purpose, NULL is returned. The reference count of the returned KSP is
223: not incremented: it should not be destroyed by the user.
225: Input Parameters:
226: . ksp - the Krylov space context
228: Output Parameters:
229: . kspest the eigenvalue estimation Krylov space context
231: Level: intermediate
233: .seealso: KSPChebyshevEstEigSet()
234: @*/
235: PetscErrorCode KSPChebyshevEstEigGetKSP(KSP ksp, KSP *kspest)
236: {
241: *kspest = NULL;
242: PetscTryMethod(ksp,"KSPChebyshevEstEigGetKSP_C",(KSP,KSP*),(ksp,kspest));
243: return(0);
244: }
248: static PetscErrorCode KSPChebyshevEstEigGetKSP_Chebyshev(KSP ksp, KSP *kspest)
249: {
250: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
253: *kspest = cheb->kspest;
254: return(0);
255: }
259: static PetscErrorCode KSPSetFromOptions_Chebyshev(PetscOptions *PetscOptionsObject,KSP ksp)
260: {
261: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
263: PetscInt neigarg = 2, nestarg = 4;
264: PetscReal eminmax[2] = {0., 0.};
265: PetscReal tform[4] = {PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE};
266: PetscBool flgeig, flgest;
269: PetscOptionsHead(PetscOptionsObject,"KSP Chebyshev Options");
270: PetscOptionsInt("-ksp_chebyshev_esteig_steps","Number of est steps in Chebyshev","",cheb->eststeps,&cheb->eststeps,NULL);
271: PetscOptionsRealArray("-ksp_chebyshev_eigenvalues","extreme eigenvalues","KSPChebyshevSetEigenvalues",eminmax,&neigarg,&flgeig);
272: if (flgeig) {
273: if (neigarg != 2) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"-ksp_chebyshev_eigenvalues: must specify 2 parameters, min and max eigenvalues");
274: KSPChebyshevSetEigenvalues(ksp, eminmax[1], eminmax[0]);
275: }
276: PetscOptionsRealArray("-ksp_chebyshev_esteig","estimate eigenvalues using a Krylov method, then use this transform for Chebyshev eigenvalue bounds","KSPChebyshevEstEigSet",tform,&nestarg,&flgest);
277: if (flgest) {
278: switch (nestarg) {
279: case 0:
280: KSPChebyshevEstEigSet(ksp,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE);
281: break;
282: case 2: /* Base everything on the max eigenvalues */
283: KSPChebyshevEstEigSet(ksp,PETSC_DECIDE,tform[0],PETSC_DECIDE,tform[1]);
284: break;
285: case 4: /* Use the full 2x2 linear transformation */
286: KSPChebyshevEstEigSet(ksp,tform[0],tform[1],tform[2],tform[3]);
287: break;
288: default: SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"Must specify either 0, 2, or 4 parameters for eigenvalue estimation");
289: }
290: }
292: /* We need to estimate eigenvalues; need to set this here so that KSPSetFromOptions() is called on the estimator */
293: if ((cheb->emin == 0. || cheb->emax == 0.) && !cheb->kspest) {
294: KSPChebyshevEstEigSet(ksp,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE);
295: }
297: if (cheb->kspest) {
298: PetscBool estrand = PETSC_FALSE;
299: PetscOptionsBool("-ksp_chebyshev_esteig_random","Use Random right hand side for eigenvalue estimation","KSPChebyshevEstEigSetRandom",estrand,&estrand,NULL);
300: if (estrand) {
301: PetscRandom random;
302: PetscRandomCreate(PetscObjectComm((PetscObject)ksp),&random);
303: PetscObjectSetOptionsPrefix((PetscObject)random,((PetscObject)ksp)->prefix);
304: PetscObjectAppendOptionsPrefix((PetscObject)random,"ksp_chebyshev_esteig_");
305: PetscRandomSetFromOptions(random);
306: KSPChebyshevEstEigSetRandom(ksp,random);
307: PetscRandomDestroy(&random);
308: }
309: }
311: if (cheb->kspest) {
312: KSPSetFromOptions(cheb->kspest);
313: }
314: PetscOptionsTail();
315: return(0);
316: }
320: /*
321: * Must be passed a KSP solver that has "converged", with KSPSetComputeEigenvalues() called before the solve
322: */
323: static PetscErrorCode KSPChebyshevComputeExtremeEigenvalues_Private(KSP kspest,PetscReal *emin,PetscReal *emax)
324: {
326: PetscInt n,neig;
327: PetscReal *re,*im,min,max;
330: KSPGetIterationNumber(kspest,&n);
331: PetscMalloc2(n,&re,n,&im);
332: KSPComputeEigenvalues(kspest,n,re,im,&neig);
333: min = PETSC_MAX_REAL;
334: max = PETSC_MIN_REAL;
335: for (n=0; n<neig; n++) {
336: min = PetscMin(min,re[n]);
337: max = PetscMax(max,re[n]);
338: }
339: PetscFree2(re,im);
340: *emax = max;
341: *emin = min;
342: return(0);
343: }
347: static PetscErrorCode KSPSolve_Chebyshev(KSP ksp)
348: {
349: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
351: PetscInt k,kp1,km1,maxit,ktmp,i;
352: PetscScalar alpha,omegaprod,mu,omega,Gamma,c[3],scale;
353: PetscReal rnorm = 0.0;
354: Vec sol_orig,b,p[3],r;
355: Mat Amat,Pmat;
356: PetscBool diagonalscale;
359: PCGetDiagonalScale(ksp->pc,&diagonalscale);
360: if (diagonalscale) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Krylov method %s does not support diagonal scaling",((PetscObject)ksp)->type_name);
362: PCGetOperators(ksp->pc,&Amat,&Pmat);
363: if (cheb->kspest) {
364: PetscObjectId amatid, pmatid;
365: PetscObjectState amatstate, pmatstate;
367: PetscObjectGetId((PetscObject)Amat,&amatid);
368: PetscObjectGetId((PetscObject)Pmat,&pmatid);
369: PetscObjectStateGet((PetscObject)Amat,&amatstate);
370: PetscObjectStateGet((PetscObject)Pmat,&pmatstate);
371: if (amatid != cheb->amatid || pmatid != cheb->pmatid || amatstate != cheb->amatstate || pmatstate != cheb->pmatstate) {
372: PetscReal max=0.0,min=0.0;
373: Vec X,B;
374: X = ksp->work[0];
375: if (cheb->random) {
376: B = ksp->work[1];
377: VecSetRandom(B,cheb->random);
378: } else {
379: B = ksp->vec_rhs;
380: }
381: KSPSolve(cheb->kspest,B,X);
383: if (ksp->guess_zero) {
384: VecZeroEntries(X);
385: }
386: KSPChebyshevComputeExtremeEigenvalues_Private(cheb->kspest,&min,&max);
388: cheb->emin = cheb->tform[0]*min + cheb->tform[1]*max;
389: cheb->emax = cheb->tform[2]*min + cheb->tform[3]*max;
391: cheb->amatid = amatid;
392: cheb->pmatid = pmatid;
393: cheb->amatstate = amatstate;
394: cheb->pmatstate = pmatstate;
395: }
396: }
398: ksp->its = 0;
399: maxit = ksp->max_it;
401: /* These three point to the three active solutions, we
402: rotate these three at each solution update */
403: km1 = 0; k = 1; kp1 = 2;
404: sol_orig = ksp->vec_sol; /* ksp->vec_sol will be asigned to rotating vector p[k], thus save its address */
405: b = ksp->vec_rhs;
406: p[km1] = sol_orig;
407: p[k] = ksp->work[0];
408: p[kp1] = ksp->work[1];
409: r = ksp->work[2];
411: /* use scale*B as our preconditioner */
412: scale = 2.0/(cheb->emax + cheb->emin);
414: /* -alpha <= scale*lambda(B^{-1}A) <= alpha */
415: alpha = 1.0 - scale*(cheb->emin);
416: Gamma = 1.0;
417: mu = 1.0/alpha;
418: omegaprod = 2.0/alpha;
420: c[km1] = 1.0;
421: c[k] = mu;
423: if (!ksp->guess_zero) {
424: KSP_MatMult(ksp,Amat,p[km1],r); /* r = b - A*p[km1] */
425: VecAYPX(r,-1.0,b);
426: } else {
427: VecCopy(b,r);
428: }
430: KSP_PCApply(ksp,r,p[k]); /* p[k] = scale B^{-1}r + p[km1] */
431: VecAYPX(p[k],scale,p[km1]);
433: for (i=0; i<maxit; i++) {
434: PetscObjectSAWsTakeAccess((PetscObject)ksp);
436: ksp->its++;
437: PetscObjectSAWsGrantAccess((PetscObject)ksp);
438: c[kp1] = 2.0*mu*c[k] - c[km1];
439: omega = omegaprod*c[k]/c[kp1];
441: KSP_MatMult(ksp,Amat,p[k],r); /* r = b - Ap[k] */
442: VecAYPX(r,-1.0,b);
443: KSP_PCApply(ksp,r,p[kp1]); /* p[kp1] = B^{-1}r */
444: ksp->vec_sol = p[k];
446: /* calculate residual norm if requested */
447: if (ksp->normtype != KSP_NORM_NONE || ksp->numbermonitors) {
448: if (ksp->normtype == KSP_NORM_UNPRECONDITIONED) {
449: VecNorm(r,NORM_2,&rnorm);
450: } else {
451: VecNorm(p[kp1],NORM_2,&rnorm);
452: }
453: PetscObjectSAWsTakeAccess((PetscObject)ksp);
454: ksp->rnorm = rnorm;
455: PetscObjectSAWsGrantAccess((PetscObject)ksp);
456: KSPLogResidualHistory(ksp,rnorm);
457: KSPMonitor(ksp,i,rnorm);
458: (*ksp->converged)(ksp,i,rnorm,&ksp->reason,ksp->cnvP);
459: if (ksp->reason) break;
460: }
462: /* y^{k+1} = omega(y^{k} - y^{k-1} + Gamma*r^{k}) + y^{k-1} */
463: VecAXPBYPCZ(p[kp1],1.0-omega,omega,omega*Gamma*scale,p[km1],p[k]);
465: ktmp = km1;
466: km1 = k;
467: k = kp1;
468: kp1 = ktmp;
469: }
470: if (!ksp->reason) {
471: if (ksp->normtype != KSP_NORM_NONE) {
472: KSP_MatMult(ksp,Amat,p[k],r); /* r = b - Ap[k] */
473: VecAYPX(r,-1.0,b);
474: if (ksp->normtype == KSP_NORM_UNPRECONDITIONED) {
475: VecNorm(r,NORM_2,&rnorm);
476: } else {
477: KSP_PCApply(ksp,r,p[kp1]); /* p[kp1] = B^{-1}r */
478: VecNorm(p[kp1],NORM_2,&rnorm);
479: }
480: PetscObjectSAWsTakeAccess((PetscObject)ksp);
481: ksp->rnorm = rnorm;
482: PetscObjectSAWsGrantAccess((PetscObject)ksp);
483: ksp->vec_sol = p[k];
484: KSPLogResidualHistory(ksp,rnorm);
485: KSPMonitor(ksp,i,rnorm);
486: }
487: if (ksp->its >= ksp->max_it) {
488: if (ksp->normtype != KSP_NORM_NONE) {
489: (*ksp->converged)(ksp,i,rnorm,&ksp->reason,ksp->cnvP);
490: if (!ksp->reason) ksp->reason = KSP_DIVERGED_ITS;
491: } else ksp->reason = KSP_CONVERGED_ITS;
492: }
493: }
495: /* make sure solution is in vector x */
496: ksp->vec_sol = sol_orig;
497: if (k) {
498: VecCopy(p[k],sol_orig);
499: }
500: return(0);
501: }
505: static PetscErrorCode KSPView_Chebyshev(KSP ksp,PetscViewer viewer)
506: {
507: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
509: PetscBool iascii;
512: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
513: if (iascii) {
514: PetscViewerASCIIPrintf(viewer," Chebyshev: eigenvalue estimates: min = %g, max = %g\n",(double)cheb->emin,(double)cheb->emax);
515: if (cheb->kspest) {
516: PetscViewerASCIIPrintf(viewer," Chebyshev: eigenvalues estimated using %s with translations [%g %g; %g %g]\n",((PetscObject) cheb->kspest)->type_name,(double)cheb->tform[0],(double)cheb->tform[1],(double)cheb->tform[2],(double)cheb->tform[3]);
517: PetscViewerASCIIPushTab(viewer);
518: KSPView(cheb->kspest,viewer);
519: PetscViewerASCIIPopTab(viewer);
520: if (cheb->random) {
521: PetscViewerASCIIPrintf(viewer," Chebyshev: estimating eigenvalues using random right hand side\n");
522: PetscViewerASCIIPushTab(viewer);
523: PetscRandomView(cheb->random,viewer);
524: PetscViewerASCIIPopTab(viewer);
525: }
526: }
527: }
528: return(0);
529: }
533: static PetscErrorCode KSPDestroy_Chebyshev(KSP ksp)
534: {
535: KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;
539: KSPDestroy(&cheb->kspest);
540: PetscRandomDestroy(&cheb->random);
541: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetEigenvalues_C",NULL);
542: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevEstEigSet_C",NULL);
543: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevEstEigSetRandom_C",NULL);
544: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevEstEigGetKSP_C",NULL);
545: KSPDestroyDefault(ksp);
546: return(0);
547: }
549: /*MC
550: KSPCHEBYSHEV - The preconditioned Chebyshev iterative method
552: Options Database Keys:
553: + -ksp_chebyshev_eigenvalues <emin,emax> - set approximations to the smallest and largest eigenvalues
554: of the preconditioned operator. If these are accurate you will get much faster convergence.
555: . -ksp_chebyshev_esteig <a,b,c,d> - estimate eigenvalues using a Krylov method, then use this
556: transform for Chebyshev eigenvalue bounds (KSPChebyshevEstEigSet())
557: . -ksp_chebyshev_esteig_steps - number of estimation steps
558: + -ksp_chebyshev_esteig_random - use random right hand side for eigenvalue estimation (KSPChebyshevEstEigSetRandom())
561: Level: beginner
563: Notes: The Chebyshev method requires both the matrix and preconditioner to
564: be symmetric positive (semi) definite.
565: Only support for left preconditioning.
567: Chebyshev is configured as a smoother by default, targetting the "upper" part of the spectrum.
568: The user should call KSPChebyshevSetEigenvalues() if they have eigenvalue estimates.
570: .seealso: KSPCreate(), KSPSetType(), KSPType (for list of available types), KSP,
571: KSPChebyshevSetEigenvalues(), KSPChebyshevEstEigSet(), KSPChebyshevEstEigSetRandom(), KSPRICHARDSON, KSPCG, PCMG
573: M*/
577: PETSC_EXTERN PetscErrorCode KSPCreate_Chebyshev(KSP ksp)
578: {
580: KSP_Chebyshev *chebyshevP;
583: PetscNewLog(ksp,&chebyshevP);
585: ksp->data = (void*)chebyshevP;
586: KSPSetSupportedNorm(ksp,KSP_NORM_PRECONDITIONED,PC_LEFT,3);
587: KSPSetSupportedNorm(ksp,KSP_NORM_UNPRECONDITIONED,PC_LEFT,2);
589: chebyshevP->emin = 0.;
590: chebyshevP->emax = 0.;
592: chebyshevP->tform[0] = 0.0;
593: chebyshevP->tform[1] = 0.1;
594: chebyshevP->tform[2] = 0;
595: chebyshevP->tform[3] = 1.1;
596: chebyshevP->eststeps = 10;
597:
598: ksp->ops->setup = KSPSetUp_Chebyshev;
599: ksp->ops->solve = KSPSolve_Chebyshev;
600: ksp->ops->destroy = KSPDestroy_Chebyshev;
601: ksp->ops->buildsolution = KSPBuildSolutionDefault;
602: ksp->ops->buildresidual = KSPBuildResidualDefault;
603: ksp->ops->setfromoptions = KSPSetFromOptions_Chebyshev;
604: ksp->ops->view = KSPView_Chebyshev;
605: ksp->ops->reset = KSPReset_Chebyshev;
607: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetEigenvalues_C",KSPChebyshevSetEigenvalues_Chebyshev);
608: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevEstEigSet_C",KSPChebyshevEstEigSet_Chebyshev);
609: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevEstEigSetRandom_C",KSPChebyshevEstEigSetRandom_Chebyshev);
610: PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevEstEigGetKSP_C",KSPChebyshevEstEigGetKSP_Chebyshev);
611: return(0);
612: }