Actual source code: cheby.c

petsc-master 2015-03-03
Report Typos and Errors
  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:   }
 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:   KSPChebyshevSetEstimateEigenvalues(ksp,0.,0.,0.,0.); /* Destroy any estimation setup */
 46:   return(0);
 47: }

 51: static PetscErrorCode KSPChebyshevSetEstimateEigenvalues_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,"est_");
 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:    KSPChebyshevSetEstimateEigenvalues - 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_estimate_eigenvalues 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 KSPChebyshevSetEstimateEigenvalues(KSP ksp,PetscReal a,PetscReal b,PetscReal c,PetscReal d)
176: {

185:   PetscTryMethod(ksp,"KSPChebyshevSetEstimateEigenvalues_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:    Level: intermediate

202: .seealso: KSPChebyshevSetEstimateEigenvalues(), PetscRandomCreate()
203: @*/
204: PetscErrorCode KSPChebyshevEstEigSetRandom(KSP ksp,PetscRandom random)
205: {

211:   PetscTryMethod(ksp,"KSPChebyshevEstEigSetRandom_C",(KSP,PetscRandom),(ksp,random));
212:   return(0);
213: }

217: /*@
218:   KSPChebyshevGetEstEigKSP - Get the Krylov method context used to estimate eigenvalues for the Chebyshev method.  If
219:   a Krylov method is not being used for this purpose, NULL is returned.  The reference count of the returned KSP is
220:   not incremented: it should not be destroyed by the user.

222:   Input Parameters:
223: . ksp - the Krylov space context

225:   Output Parameters:
226: . kspest the eigenvalue estimation Krylov space context

228:   Level: intermediate

230: .seealso: KSPChebyshevSetEstimateEigenvalues()
231: @*/
232: PetscErrorCode KSPChebyshevGetEstEigKSP(KSP ksp, KSP *kspest)
233: {

238:   *kspest = NULL;
239:   PetscTryMethod(ksp,"KSPChebyshevGetEstEigKSP_C",(KSP,KSP*),(ksp,kspest));
240:   return(0);
241: }

245: static PetscErrorCode KSPChebyshevGetEstEigKSP_Chebyshev(KSP ksp, KSP *kspest)
246: {
247:   KSP_Chebyshev *cheb = (KSP_Chebyshev*)ksp->data;

250:   *kspest = cheb->kspest;
251:   return(0);
252: }

256: PetscErrorCode KSPSetFromOptions_Chebyshev(PetscOptions *PetscOptionsObject,KSP ksp)
257: {
258:   KSP_Chebyshev  *cheb = (KSP_Chebyshev*)ksp->data;
260:   PetscInt       neigarg = 2, nestarg = 4;
261:   PetscReal      eminmax[2] = {0., 0.};
262:   PetscReal      tform[4] = {PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE};
263:   PetscBool      flgeig, flgest;

266:   PetscOptionsHead(PetscOptionsObject,"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",eminmax,&neigarg,&flgeig);
269:   if (flgeig) {
270:     if (neigarg != 2) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"-ksp_chebyshev_eigenvalues: must specify 2 parameters, min and max eigenvalues");
271:     KSPChebyshevSetEigenvalues(ksp, eminmax[1], eminmax[0]);
272:   }
273:   PetscOptionsRealArray("-ksp_chebyshev_estimate_eigenvalues","estimate eigenvalues using a Krylov method, then use this transform for Chebyshev eigenvalue bounds","KSPChebyshevSetEstimateEigenvalues",tform,&nestarg,&flgest);
274:   if (flgest) {
275:     switch (nestarg) {
276:     case 0:
277:       KSPChebyshevSetEstimateEigenvalues(ksp,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE);
278:       break;
279:     case 2:                     /* Base everything on the max eigenvalues */
280:       KSPChebyshevSetEstimateEigenvalues(ksp,PETSC_DECIDE,tform[0],PETSC_DECIDE,tform[1]);
281:       break;
282:     case 4:                     /* Use the full 2x2 linear transformation */
283:       KSPChebyshevSetEstimateEigenvalues(ksp,tform[0],tform[1],tform[2],tform[3]);
284:       break;
285:     default: SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_INCOMP,"Must specify either 0, 2, or 4 parameters for eigenvalue estimation");
286:     }
287:   }

289:   if (cheb->kspest) {
290:     PetscBool estrand = PETSC_FALSE;
291:     PetscOptionsBool("-ksp_chebyshev_estimate_eigenvalues_random","Use Random right hand side for eigenvalue estimation","KSPChebyshevEstEigSetRandom",estrand,&estrand,NULL);
292:     if (estrand) {
293:       PetscRandom random;
294:       PetscRandomCreate(PetscObjectComm((PetscObject)ksp),&random);
295:       PetscObjectSetOptionsPrefix((PetscObject)random,((PetscObject)ksp)->prefix);
296:       PetscObjectAppendOptionsPrefix((PetscObject)random,"ksp_chebyshev_estimate_eigenvalues_");
297:       PetscRandomSetFromOptions(random);
298:       KSPChebyshevEstEigSetRandom(ksp,random);
299:       PetscRandomDestroy(&random);
300:     }
301:   }

303:   if (cheb->kspest) {
304:     KSPSetFromOptions(cheb->kspest);
305:   }
306:   PetscOptionsTail();
307:   return(0);
308: }

312: /*
313:  * Must be passed a KSP solver that has "converged", with KSPSetComputeEigenvalues() called before the solve
314:  */
315: static PetscErrorCode KSPChebyshevComputeExtremeEigenvalues_Private(KSP kspest,PetscReal *emin,PetscReal *emax)
316: {
318:   PetscInt       n,neig;
319:   PetscReal      *re,*im,min,max;

322:   KSPGetIterationNumber(kspest,&n);
323:   PetscMalloc2(n,&re,n,&im);
324:   KSPComputeEigenvalues(kspest,n,re,im,&neig);
325:   min  = PETSC_MAX_REAL;
326:   max  = PETSC_MIN_REAL;
327:   for (n=0; n<neig; n++) {
328:     min = PetscMin(min,re[n]);
329:     max = PetscMax(max,re[n]);
330:   }
331:   PetscFree2(re,im);
332:   *emax = max;
333:   *emin = min;
334:   return(0);
335: }

339: PetscErrorCode KSPSolve_Chebyshev(KSP ksp)
340: {
341:   KSP_Chebyshev  *cheb = (KSP_Chebyshev*)ksp->data;
343:   PetscInt       k,kp1,km1,maxit,ktmp,i;
344:   PetscScalar    alpha,omegaprod,mu,omega,Gamma,c[3],scale;
345:   PetscReal      rnorm = 0.0;
346:   Vec            sol_orig,b,p[3],r;
347:   Mat            Amat,Pmat;
348:   PetscBool      diagonalscale;

351:   PCGetDiagonalScale(ksp->pc,&diagonalscale);
352:   if (diagonalscale) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"Krylov method %s does not support diagonal scaling",((PetscObject)ksp)->type_name);

354:   PCGetOperators(ksp->pc,&Amat,&Pmat);
355:   if (cheb->kspest) {
356:     PetscObjectId    amatid,    pmatid;
357:     PetscObjectState amatstate, pmatstate;

359:     PetscObjectGetId((PetscObject)Amat,&amatid);
360:     PetscObjectGetId((PetscObject)Pmat,&pmatid);
361:     PetscObjectStateGet((PetscObject)Amat,&amatstate);
362:     PetscObjectStateGet((PetscObject)Pmat,&pmatstate);
363:     if (amatid != cheb->amatid || pmatid != cheb->pmatid || amatstate != cheb->amatstate || pmatstate != cheb->pmatstate) {
364:       PetscReal max=0.0,min=0.0;
365:       Vec       X,B;
366:       X = ksp->work[0];
367:       if (cheb->random) {
368:         B    = ksp->work[1];
369:         VecSetRandom(B,cheb->random);
370:       } else {
371:         B = ksp->vec_rhs;
372:       }
373:       KSPSolve(cheb->kspest,B,X);

375:       if (ksp->guess_zero) {
376:         VecZeroEntries(X);
377:       }
378:       KSPChebyshevComputeExtremeEigenvalues_Private(cheb->kspest,&min,&max);

380:       cheb->emin = cheb->tform[0]*min + cheb->tform[1]*max;
381:       cheb->emax = cheb->tform[2]*min + cheb->tform[3]*max;

383:       cheb->amatid    = amatid;
384:       cheb->pmatid    = pmatid;
385:       cheb->amatstate = amatstate;
386:       cheb->pmatstate = pmatstate;
387:     }
388:   }

390:   ksp->its = 0;
391:   maxit    = ksp->max_it;

393:   /* These three point to the three active solutions, we
394:      rotate these three at each solution update */
395:   km1      = 0; k = 1; kp1 = 2;
396:   sol_orig = ksp->vec_sol; /* ksp->vec_sol will be asigned to rotating vector p[k], thus save its address */
397:   b        = ksp->vec_rhs;
398:   p[km1]   = sol_orig;
399:   p[k]     = ksp->work[0];
400:   p[kp1]   = ksp->work[1];
401:   r        = ksp->work[2];

403:   /* use scale*B as our preconditioner */
404:   scale = 2.0/(cheb->emax + cheb->emin);

406:   /*   -alpha <=  scale*lambda(B^{-1}A) <= alpha   */
407:   alpha     = 1.0 - scale*(cheb->emin);
408:   Gamma     = 1.0;
409:   mu        = 1.0/alpha;
410:   omegaprod = 2.0/alpha;

412:   c[km1] = 1.0;
413:   c[k]   = mu;

415:   if (!ksp->guess_zero) {
416:     KSP_MatMult(ksp,Amat,p[km1],r);     /*  r = b - A*p[km1] */
417:     VecAYPX(r,-1.0,b);
418:   } else {
419:     VecCopy(b,r);
420:   }

422:   KSP_PCApply(ksp,r,p[k]);  /* p[k] = scale B^{-1}r + p[km1] */
423:   VecAYPX(p[k],scale,p[km1]);

425:   for (i=0; i<maxit; i++) {
426:     PetscObjectSAWsTakeAccess((PetscObject)ksp);

428:     ksp->its++;
429:     PetscObjectSAWsGrantAccess((PetscObject)ksp);
430:     c[kp1] = 2.0*mu*c[k] - c[km1];
431:     omega  = omegaprod*c[k]/c[kp1];

433:     KSP_MatMult(ksp,Amat,p[k],r);          /*  r = b - Ap[k]    */
434:     VecAYPX(r,-1.0,b);
435:     KSP_PCApply(ksp,r,p[kp1]);             /*  p[kp1] = B^{-1}r  */
436:     ksp->vec_sol = p[k];

438:     /* calculate residual norm if requested */
439:     if (ksp->normtype != KSP_NORM_NONE || ksp->numbermonitors) {
440:       if (ksp->normtype == KSP_NORM_UNPRECONDITIONED) {
441:         VecNorm(r,NORM_2,&rnorm);
442:       } else {
443:         VecNorm(p[kp1],NORM_2,&rnorm);
444:       }
445:       PetscObjectSAWsTakeAccess((PetscObject)ksp);
446:       ksp->rnorm   = rnorm;
447:       PetscObjectSAWsGrantAccess((PetscObject)ksp);
448:       KSPLogResidualHistory(ksp,rnorm);
449:       KSPMonitor(ksp,i,rnorm);
450:       (*ksp->converged)(ksp,i,rnorm,&ksp->reason,ksp->cnvP);
451:       if (ksp->reason) break;
452:     }

454:     /* y^{k+1} = omega(y^{k} - y^{k-1} + Gamma*r^{k}) + y^{k-1} */
455:     VecAXPBYPCZ(p[kp1],1.0-omega,omega,omega*Gamma*scale,p[km1],p[k]);

457:     ktmp = km1;
458:     km1  = k;
459:     k    = kp1;
460:     kp1  = ktmp;
461:   }
462:   if (!ksp->reason) {
463:     if (ksp->normtype != KSP_NORM_NONE) {
464:       KSP_MatMult(ksp,Amat,p[k],r);       /*  r = b - Ap[k]    */
465:       VecAYPX(r,-1.0,b);
466:       if (ksp->normtype == KSP_NORM_UNPRECONDITIONED) {
467:         VecNorm(r,NORM_2,&rnorm);
468:       } else {
469:         KSP_PCApply(ksp,r,p[kp1]); /* p[kp1] = B^{-1}r */
470:         VecNorm(p[kp1],NORM_2,&rnorm);
471:       }
472:       PetscObjectSAWsTakeAccess((PetscObject)ksp);
473:       ksp->rnorm   = rnorm;
474:       PetscObjectSAWsGrantAccess((PetscObject)ksp);
475:       ksp->vec_sol = p[k];
476:       KSPLogResidualHistory(ksp,rnorm);
477:       KSPMonitor(ksp,i,rnorm);
478:     }
479:     if (ksp->its >= ksp->max_it) {
480:       if (ksp->normtype != KSP_NORM_NONE) {
481:         (*ksp->converged)(ksp,i,rnorm,&ksp->reason,ksp->cnvP);
482:         if (!ksp->reason) ksp->reason = KSP_DIVERGED_ITS;
483:       } else ksp->reason = KSP_CONVERGED_ITS;
484:     }
485:   }

487:   /* make sure solution is in vector x */
488:   ksp->vec_sol = sol_orig;
489:   if (k) {
490:     VecCopy(p[k],sol_orig);
491:   }
492:   return(0);
493: }

497: PetscErrorCode KSPView_Chebyshev(KSP ksp,PetscViewer viewer)
498: {
499:   KSP_Chebyshev  *cheb = (KSP_Chebyshev*)ksp->data;
501:   PetscBool      iascii;

504:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
505:   if (iascii) {
506:     PetscViewerASCIIPrintf(viewer,"  Chebyshev: eigenvalue estimates:  min = %g, max = %g\n",(double)cheb->emin,(double)cheb->emax);
507:     if (cheb->kspest) {
508:       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]);
509:       PetscViewerASCIIPushTab(viewer);
510:       KSPView(cheb->kspest,viewer);
511:       PetscViewerASCIIPopTab(viewer);
512:       if (cheb->random) {
513:         PetscViewerASCIIPrintf(viewer,"  Chebyshev: estimating eigenvalues using random right hand side\n");
514:         PetscViewerASCIIPushTab(viewer);
515:         PetscRandomView(cheb->random,viewer);
516:         PetscViewerASCIIPopTab(viewer);
517:       }
518:     }
519:   }
520:   return(0);
521: }

525: PetscErrorCode KSPDestroy_Chebyshev(KSP ksp)
526: {
527:   KSP_Chebyshev  *cheb = (KSP_Chebyshev*)ksp->data;

531:   KSPDestroy(&cheb->kspest);
532:   PetscRandomDestroy(&cheb->random);
533:   PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetEigenvalues_C",NULL);
534:   PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetEstimateEigenvalues_C",NULL);
535:   PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevEstEigSetRandom_C",NULL);
536:   PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevGetEstEigKSP_C",NULL);
537:   KSPDestroyDefault(ksp);
538:   return(0);
539: }

541: /*MC
542:      KSPCHEBYSHEV - The preconditioned Chebyshev iterative method

544:    Options Database Keys:
545: +   -ksp_chebyshev_eigenvalues <emin,emax> - set approximations to the smallest and largest eigenvalues
546:                   of the preconditioned operator. If these are accurate you will get much faster convergence.
547: -   -ksp_chebyshev_estimate_eigenvalues <a,b,c,d> - estimate eigenvalues using a Krylov method, then use this
548:                   transform for Chebyshev eigenvalue bounds (KSPChebyshevSetEstimateEigenvalues)


551:    Level: beginner

553:    Notes: The Chebyshev method requires both the matrix and preconditioner to
554:           be symmetric positive (semi) definite.
555:           Only support for left preconditioning.

557:           Chebyshev is configured as a smoother by default, targetting the "upper" part of the spectrum.
558:           The user should call KSPChebyshevSetEigenvalues() if they have eigenvalue estimates.

560: .seealso:  KSPCreate(), KSPSetType(), KSPType (for list of available types), KSP,
561:            KSPChebyshevSetEigenvalues(), KSPChebyshevSetEstimateEigenvalues(), KSPRICHARDSON, KSPCG, PCMG

563: M*/

567: PETSC_EXTERN PetscErrorCode KSPCreate_Chebyshev(KSP ksp)
568: {
570:   KSP_Chebyshev  *chebyshevP;

573:   PetscNewLog(ksp,&chebyshevP);

575:   ksp->data = (void*)chebyshevP;
576:   KSPSetSupportedNorm(ksp,KSP_NORM_PRECONDITIONED,PC_LEFT,3);
577:   KSPSetSupportedNorm(ksp,KSP_NORM_UNPRECONDITIONED,PC_LEFT,2);

579:   chebyshevP->emin = 0.;
580:   chebyshevP->emax = 0.;

582:   chebyshevP->tform[0] = 0.0;
583:   chebyshevP->tform[1] = 0.1;
584:   chebyshevP->tform[2] = 0;
585:   chebyshevP->tform[3] = 1.1;
586:   chebyshevP->eststeps = 10;
587: 
588:   ksp->ops->setup          = KSPSetUp_Chebyshev;
589:   ksp->ops->solve          = KSPSolve_Chebyshev;
590:   ksp->ops->destroy        = KSPDestroy_Chebyshev;
591:   ksp->ops->buildsolution  = KSPBuildSolutionDefault;
592:   ksp->ops->buildresidual  = KSPBuildResidualDefault;
593:   ksp->ops->setfromoptions = KSPSetFromOptions_Chebyshev;
594:   ksp->ops->view           = KSPView_Chebyshev;
595:   ksp->ops->reset          = KSPReset_Chebyshev;

597:   PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetEigenvalues_C",KSPChebyshevSetEigenvalues_Chebyshev);
598:   PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevSetEstimateEigenvalues_C",KSPChebyshevSetEstimateEigenvalues_Chebyshev);
599:   PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevEstEigSetRandom_C",KSPChebyshevEstEigSetRandom_Chebyshev);
600:   PetscObjectComposeFunction((PetscObject)ksp,"KSPChebyshevGetEstEigKSP_C",KSPChebyshevGetEstEigKSP_Chebyshev);
601:   return(0);
602: }