Actual source code: fcg.c

  1: /*
  2:     This file implements the FCG (Flexible Conjugate Gradient) method

  4: */

  6: #include <../src/ksp/ksp/impls/fcg/fcgimpl.h>
  7: extern PetscErrorCode KSPComputeExtremeSingularValues_CG(KSP,PetscReal*,PetscReal*);
  8: extern PetscErrorCode KSPComputeEigenvalues_CG(KSP,PetscInt,PetscReal*,PetscReal*,PetscInt*);

 10: #define KSPFCG_DEFAULT_MMAX 30          /* maximum number of search directions to keep */
 11: #define KSPFCG_DEFAULT_NPREALLOC 10     /* number of search directions to preallocate */
 12: #define KSPFCG_DEFAULT_VECB 5           /* number of search directions to allocate each time new direction vectors are needed */
 13: #define KSPFCG_DEFAULT_TRUNCSTRAT KSP_FCD_TRUNC_TYPE_NOTAY

 15: static PetscErrorCode KSPAllocateVectors_FCG(KSP ksp, PetscInt nvecsneeded, PetscInt chunksize)
 16: {
 17:   PetscInt        i;
 18:   KSP_FCG         *fcg = (KSP_FCG*)ksp->data;
 19:   PetscInt        nnewvecs, nvecsprev;

 21:   /* Allocate enough new vectors to add chunksize new vectors, reach nvecsneedtotal, or to reach mmax+1, whichever is smallest */
 22:   if (fcg->nvecs < PetscMin(fcg->mmax+1,nvecsneeded)) {
 23:     nvecsprev = fcg->nvecs;
 24:     nnewvecs = PetscMin(PetscMax(nvecsneeded-fcg->nvecs,chunksize),fcg->mmax+1-fcg->nvecs);
 25:     KSPCreateVecs(ksp,nnewvecs,&fcg->pCvecs[fcg->nchunks],0,NULL);
 26:     PetscLogObjectParents((PetscObject)ksp,nnewvecs,fcg->pCvecs[fcg->nchunks]);
 27:     KSPCreateVecs(ksp,nnewvecs,&fcg->pPvecs[fcg->nchunks],0,NULL);
 28:     PetscLogObjectParents((PetscObject)ksp,nnewvecs,fcg->pPvecs[fcg->nchunks]);
 29:     fcg->nvecs += nnewvecs;
 30:     for (i=0;i<nnewvecs;++i) {
 31:       fcg->Cvecs[nvecsprev + i] = fcg->pCvecs[fcg->nchunks][i];
 32:       fcg->Pvecs[nvecsprev + i] = fcg->pPvecs[fcg->nchunks][i];
 33:     }
 34:     fcg->chunksizes[fcg->nchunks] = nnewvecs;
 35:     ++fcg->nchunks;
 36:   }
 37:   return 0;
 38: }

 40: static PetscErrorCode    KSPSetUp_FCG(KSP ksp)
 41: {
 42:   KSP_FCG        *fcg = (KSP_FCG*)ksp->data;
 43:   PetscInt       maxit = ksp->max_it;
 44:   const PetscInt nworkstd = 2;


 47:   /* Allocate "standard" work vectors (not including the basis and transformed basis vectors) */
 48:   KSPSetWorkVecs(ksp,nworkstd);

 50:   /* Allocated space for pointers to additional work vectors
 51:    note that mmax is the number of previous directions, so we add 1 for the current direction,
 52:    and an extra 1 for the prealloc (which might be empty) */
 53:   PetscMalloc5(fcg->mmax+1,&fcg->Pvecs,fcg->mmax+1,&fcg->Cvecs,fcg->mmax+1,&fcg->pPvecs,fcg->mmax+1,&fcg->pCvecs,fcg->mmax+2,&fcg->chunksizes);
 54:   PetscLogObjectMemory((PetscObject)ksp,2*(fcg->mmax+1)*sizeof(Vec*) + 2*(fcg->mmax + 1)*sizeof(Vec**) + (fcg->mmax + 2)*sizeof(PetscInt));

 56:   /* If the requested number of preallocated vectors is greater than mmax reduce nprealloc */
 57:   if (fcg->nprealloc > fcg->mmax+1) {
 58:     PetscInfo(NULL,"Requested nprealloc=%d is greater than m_max+1=%d. Resetting nprealloc = m_max+1.\n",fcg->nprealloc, fcg->mmax+1);
 59:   }

 61:   /* Preallocate additional work vectors */
 62:   KSPAllocateVectors_FCG(ksp,fcg->nprealloc,fcg->nprealloc);
 63:   /*
 64:   If user requested computations of eigenvalues then allocate work
 65:   work space needed
 66:   */
 67:   if (ksp->calc_sings) {
 68:     /* get space to store tridiagonal matrix for Lanczos */
 69:     PetscMalloc4(maxit,&fcg->e,maxit,&fcg->d,maxit,&fcg->ee,maxit,&fcg->dd);
 70:     PetscLogObjectMemory((PetscObject)ksp,2*(maxit+1)*(sizeof(PetscScalar)+sizeof(PetscReal)));

 72:     ksp->ops->computeextremesingularvalues = KSPComputeExtremeSingularValues_CG;
 73:     ksp->ops->computeeigenvalues           = KSPComputeEigenvalues_CG;
 74:   }
 75:   return 0;
 76: }

 78: static PetscErrorCode KSPSolve_FCG(KSP ksp)
 79: {
 80:   PetscInt       i,k,idx,mi;
 81:   KSP_FCG        *fcg = (KSP_FCG*)ksp->data;
 82:   PetscScalar    alpha=0.0,beta = 0.0,dpi,s;
 83:   PetscReal      dp=0.0;
 84:   Vec            B,R,Z,X,Pcurr,Ccurr;
 85:   Mat            Amat,Pmat;
 86:   PetscInt       eigs = ksp->calc_sings; /* Variables for eigen estimation - START*/
 87:   PetscInt       stored_max_it = ksp->max_it;
 88:   PetscScalar    alphaold = 0,betaold = 1.0,*e = NULL,*d = NULL;/* Variables for eigen estimation  - FINISH */


 91: #define VecXDot(x,y,a) (((fcg->type) == (KSP_CG_HERMITIAN)) ? VecDot(x,y,a) : VecTDot(x,y,a))
 92: #define VecXMDot(a,b,c,d) (((fcg->type) == (KSP_CG_HERMITIAN)) ? VecMDot(a,b,c,d) : VecMTDot(a,b,c,d))

 94:   X             = ksp->vec_sol;
 95:   B             = ksp->vec_rhs;
 96:   R             = ksp->work[0];
 97:   Z             = ksp->work[1];

 99:   PCGetOperators(ksp->pc,&Amat,&Pmat);
100:   if (eigs) {e = fcg->e; d = fcg->d; e[0] = 0.0; }
101:   /* Compute initial residual needed for convergence check*/
102:   ksp->its = 0;
103:   if (!ksp->guess_zero) {
104:     KSP_MatMult(ksp,Amat,X,R);
105:     VecAYPX(R,-1.0,B);                    /*   r <- b - Ax     */
106:   } else {
107:     VecCopy(B,R);                         /*   r <- b (x is 0) */
108:   }
109:   switch (ksp->normtype) {
110:     case KSP_NORM_PRECONDITIONED:
111:       KSP_PCApply(ksp,R,Z);               /*   z <- Br         */
112:       VecNorm(Z,NORM_2,&dp);              /*   dp <- dqrt(z'*z) = sqrt(e'*A'*B'*B*A*e)     */
113:       KSPCheckNorm(ksp,dp);
114:       break;
115:     case KSP_NORM_UNPRECONDITIONED:
116:       VecNorm(R,NORM_2,&dp);              /*   dp <- sqrt(r'*r) = sqrt(e'*A'*A*e)     */
117:       KSPCheckNorm(ksp,dp);
118:       break;
119:     case KSP_NORM_NATURAL:
120:       KSP_PCApply(ksp,R,Z);               /*   z <- Br         */
121:       VecXDot(R,Z,&s);
122:       KSPCheckDot(ksp,s);
123:       dp = PetscSqrtReal(PetscAbsScalar(s));                   /*   dp <- sqrt(r'*z) = sqrt(e'*A'*B*A*e)  */
124:       break;
125:     case KSP_NORM_NONE:
126:       dp = 0.0;
127:       break;
128:     default: SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"%s",KSPNormTypes[ksp->normtype]);
129:   }

131:   /* Initial Convergence Check */
132:   KSPLogResidualHistory(ksp,dp);
133:   KSPMonitor(ksp,0,dp);
134:   ksp->rnorm = dp;
135:   if (ksp->normtype == KSP_NORM_NONE) {
136:     KSPConvergedSkip(ksp,0,dp,&ksp->reason,ksp->cnvP);
137:   } else {
138:     (*ksp->converged)(ksp,0,dp,&ksp->reason,ksp->cnvP);
139:   }
140:   if (ksp->reason) return 0;

142:   /* Apply PC if not already done for convergence check */
143:   if (ksp->normtype == KSP_NORM_UNPRECONDITIONED || ksp->normtype == KSP_NORM_NONE) {
144:     KSP_PCApply(ksp,R,Z);               /*   z <- Br         */
145:   }

147:   i = 0;
148:   do {
149:     ksp->its = i+1;

151:     /*  If needbe, allocate a new chunk of vectors in P and C */
152:     KSPAllocateVectors_FCG(ksp,i+1,fcg->vecb);

154:     /* Note that we wrap around and start clobbering old vectors */
155:     idx = i % (fcg->mmax+1);
156:     Pcurr = fcg->Pvecs[idx];
157:     Ccurr = fcg->Cvecs[idx];

159:     /* number of old directions to orthogonalize against */
160:     switch(fcg->truncstrat) {
161:       case KSP_FCD_TRUNC_TYPE_STANDARD:
162:         mi = fcg->mmax;
163:         break;
164:       case KSP_FCD_TRUNC_TYPE_NOTAY:
165:         mi = ((i-1) % fcg->mmax)+1;
166:         break;
167:       default:
168:         SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Unrecognized Truncation Strategy");
169:     }

171:     /* Compute a new column of P (Currently does not support modified G-S or iterative refinement)*/
172:     VecCopy(Z,Pcurr);

174:     {
175:       PetscInt l,ndots;

177:       l = PetscMax(0,i-mi);
178:       ndots = i-l;
179:       if (ndots) {
180:         PetscInt    j;
181:         Vec         *Pold,  *Cold;
182:         PetscScalar *dots;

184:         PetscMalloc3(ndots,&dots,ndots,&Cold,ndots,&Pold);
185:         for (k=l,j=0;j<ndots;++k,++j) {
186:           idx = k % (fcg->mmax+1);
187:           Cold[j] = fcg->Cvecs[idx];
188:           Pold[j] = fcg->Pvecs[idx];
189:         }
190:         VecXMDot(Z,ndots,Cold,dots);
191:         for (k=0;k<ndots;++k) {
192:           dots[k] = -dots[k];
193:         }
194:         VecMAXPY(Pcurr,ndots,dots,Pold);
195:         PetscFree3(dots,Cold,Pold);
196:       }
197:     }

199:     /* Update X and R */
200:     betaold = beta;
201:     VecXDot(Pcurr,R,&beta);                 /*  beta <- pi'*r       */
202:     KSPCheckDot(ksp,beta);
203:     KSP_MatMult(ksp,Amat,Pcurr,Ccurr);      /*  w <- A*pi (stored in ci)   */
204:     VecXDot(Pcurr,Ccurr,&dpi);              /*  dpi <- pi'*w        */
205:     alphaold = alpha;
206:     alpha = beta / dpi;                                          /*  alpha <- beta/dpi    */
207:     VecAXPY(X,alpha,Pcurr);                 /*  x <- x + alpha * pi  */
208:     VecAXPY(R,-alpha,Ccurr);                /*  r <- r - alpha * wi  */

210:     /* Compute norm for convergence check */
211:     switch (ksp->normtype) {
212:       case KSP_NORM_PRECONDITIONED:
213:         KSP_PCApply(ksp,R,Z);               /*   z <- Br             */
214:         VecNorm(Z,NORM_2,&dp);              /*   dp <- sqrt(z'*z) = sqrt(e'*A'*B'*B*A*e)  */
215:         KSPCheckNorm(ksp,dp);
216:       break;
217:       case KSP_NORM_UNPRECONDITIONED:
218:         VecNorm(R,NORM_2,&dp);              /*   dp <- sqrt(r'*r) = sqrt(e'*A'*A*e)   */
219:         KSPCheckNorm(ksp,dp);
220:         break;
221:       case KSP_NORM_NATURAL:
222:         KSP_PCApply(ksp,R,Z);               /*   z <- Br             */
223:         VecXDot(R,Z,&s);
224:         KSPCheckDot(ksp,s);
225:         dp = PetscSqrtReal(PetscAbsScalar(s));                   /*   dp <- sqrt(r'*z) = sqrt(e'*A'*B*A*e)  */
226:         break;
227:       case KSP_NORM_NONE:
228:         dp = 0.0;
229:         break;
230:       default: SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP,"%s",KSPNormTypes[ksp->normtype]);
231:     }

233:     /* Check for convergence */
234:     ksp->rnorm = dp;
235:     KSPLogResidualHistory(ksp,dp);
236:     KSPMonitor(ksp,i+1,dp);
237:     (*ksp->converged)(ksp,i+1,dp,&ksp->reason,ksp->cnvP);
238:     if (ksp->reason) break;

240:     /* Apply PC if not already done for convergence check */
241:     if (ksp->normtype == KSP_NORM_UNPRECONDITIONED || ksp->normtype == KSP_NORM_NONE) {
242:       KSP_PCApply(ksp,R,Z);               /*   z <- Br         */
243:     }

245:     /* Compute current C (which is W/dpi) */
246:     VecScale(Ccurr,1.0/dpi);              /*   w <- ci/dpi   */

248:     if (eigs) {
249:       if (i > 0) {
251:         e[i] = PetscSqrtReal(PetscAbsScalar(beta/betaold))/alphaold;
252:         d[i] = PetscSqrtReal(PetscAbsScalar(beta/betaold))*e[i] + 1.0/alpha;
253:       } else {
254:         d[i] = PetscSqrtReal(PetscAbsScalar(beta))*e[i] + 1.0/alpha;
255:       }
256:       fcg->ned = ksp->its-1;
257:     }
258:     ++i;
259:   } while (i<ksp->max_it);
260:   if (i >= ksp->max_it) ksp->reason = KSP_DIVERGED_ITS;
261:   return 0;
262: }

264: static PetscErrorCode KSPDestroy_FCG(KSP ksp)
265: {
266:   PetscInt       i;
267:   KSP_FCG        *fcg = (KSP_FCG*)ksp->data;


270:   /* Destroy "standard" work vecs */
271:   VecDestroyVecs(ksp->nwork,&ksp->work);

273:   /* Destroy P and C vectors and the arrays that manage pointers to them */
274:   if (fcg->nvecs) {
275:     for (i=0;i<fcg->nchunks;++i) {
276:       VecDestroyVecs(fcg->chunksizes[i],&fcg->pPvecs[i]);
277:       VecDestroyVecs(fcg->chunksizes[i],&fcg->pCvecs[i]);
278:     }
279:   }
280:   PetscFree5(fcg->Pvecs,fcg->Cvecs,fcg->pPvecs,fcg->pCvecs,fcg->chunksizes);
281:   /* free space used for singular value calculations */
282:   if (ksp->calc_sings) {
283:     PetscFree4(fcg->e,fcg->d,fcg->ee,fcg->dd);
284:   }
285:   KSPDestroyDefault(ksp);
286:   return 0;
287: }

289: static PetscErrorCode KSPView_FCG(KSP ksp,PetscViewer viewer)
290: {
291:   KSP_FCG        *fcg = (KSP_FCG*)ksp->data;
292:   PetscBool      iascii,isstring;
293:   const char     *truncstr;

295:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
296:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);

298:   if (fcg->truncstrat == KSP_FCD_TRUNC_TYPE_STANDARD) truncstr = "Using standard truncation strategy";
299:   else if (fcg->truncstrat == KSP_FCD_TRUNC_TYPE_NOTAY) truncstr = "Using Notay's truncation strategy";
300:   else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Undefined FCG truncation strategy");

302:   if (iascii) {
303:     PetscViewerASCIIPrintf(viewer,"  m_max=%D\n",fcg->mmax);
304:     PetscViewerASCIIPrintf(viewer,"  preallocated %D directions\n",PetscMin(fcg->nprealloc,fcg->mmax+1));
305:     PetscViewerASCIIPrintf(viewer,"  %s\n",truncstr);
306:   } else if (isstring) {
307:     PetscViewerStringSPrintf(viewer,"m_max %D nprealloc %D %s",fcg->mmax,fcg->nprealloc,truncstr);
308:   }
309:   return 0;
310: }

312: /*@
313:   KSPFCGSetMmax - set the maximum number of previous directions FCG will store for orthogonalization

315:   Note: mmax + 1 directions are stored (mmax previous ones along with a current one)
316:   and whether all are used in each iteration also depends on the truncation strategy
317:   (see KSPFCGSetTruncationType())

319:   Logically Collective on ksp

321:   Input Parameters:
322: +  ksp - the Krylov space context
323: -  mmax - the maximum number of previous directions to orthogonalize againt

325:   Level: intermediate

327: .seealso: KSPFCG, KSPFCGGetTruncationType(), KSPFCGGetNprealloc()
328: @*/
329: PetscErrorCode KSPFCGSetMmax(KSP ksp,PetscInt mmax)
330: {
331:   KSP_FCG *fcg = (KSP_FCG*)ksp->data;

335:   fcg->mmax = mmax;
336:   return 0;
337: }

339: /*@
340:   KSPFCGGetMmax - get the maximum number of previous directions FCG will store

342:   Note: FCG stores mmax+1 directions at most (mmax previous ones, and one current one)

344:    Not Collective

346:    Input Parameter:
347: .  ksp - the Krylov space context

349:    Output Parameter:
350: .  mmax - the maximum number of previous directions allowed for orthogonalization

352:    Level: intermediate

354: .seealso: KSPFCG, KSPFCGGetTruncationType(), KSPFCGGetNprealloc(), KSPFCGSetMmax()
355: @*/

357: PetscErrorCode KSPFCGGetMmax(KSP ksp,PetscInt *mmax)
358: {
359:   KSP_FCG *fcg=(KSP_FCG*)ksp->data;

362:   *mmax = fcg->mmax;
363:   return 0;
364: }

366: /*@
367:   KSPFCGSetNprealloc - set the number of directions to preallocate with FCG

369:   Logically Collective on ksp

371:   Input Parameters:
372: +  ksp - the Krylov space context
373: -  nprealloc - the number of vectors to preallocate

375:   Level: advanced

377:   Options Database:
378: . -ksp_fcg_nprealloc <N> - number of directions to preallocate

380: .seealso: KSPFCG, KSPFCGGetTruncationType(), KSPFCGGetNprealloc()
381: @*/
382: PetscErrorCode KSPFCGSetNprealloc(KSP ksp,PetscInt nprealloc)
383: {
384:   KSP_FCG *fcg=(KSP_FCG*)ksp->data;

389:   fcg->nprealloc = nprealloc;
390:   return 0;
391: }

393: /*@
394:   KSPFCGGetNprealloc - get the number of directions preallocate by FCG

396:    Not Collective

398:    Input Parameter:
399: .  ksp - the Krylov space context

401:    Output Parameter:
402: .  nprealloc - the number of directions preallocated

404:    Level: advanced

406: .seealso: KSPFCG, KSPFCGGetTruncationType(), KSPFCGSetNprealloc()
407: @*/
408: PetscErrorCode KSPFCGGetNprealloc(KSP ksp,PetscInt *nprealloc)
409: {
410:   KSP_FCG *fcg=(KSP_FCG*)ksp->data;

413:   *nprealloc = fcg->nprealloc;
414:   return 0;
415: }

417: /*@
418:   KSPFCGSetTruncationType - specify how many of its stored previous directions FCG uses during orthoganalization

420:   Logically Collective on ksp

422:   KSP_FCD_TRUNC_TYPE_STANDARD uses all (up to mmax) stored directions
423:   KSP_FCD_TRUNC_TYPE_NOTAY uses the last max(1,mod(i,mmax)) stored directions at iteration i=0,1,..

425:   Input Parameters:
426: +  ksp - the Krylov space context
427: -  truncstrat - the choice of strategy

429:   Level: intermediate

431:   Options Database:
432: . -ksp_fcg_truncation_type <standard, notay> - specify how many of its stored previous directions FCG uses during orthoganalization

434:   .seealso: KSPFCDTruncationType, KSPFCGGetTruncationType
435: @*/
436: PetscErrorCode KSPFCGSetTruncationType(KSP ksp,KSPFCDTruncationType truncstrat)
437: {
438:   KSP_FCG *fcg=(KSP_FCG*)ksp->data;

442:   fcg->truncstrat=truncstrat;
443:   return 0;
444: }

446: /*@
447:   KSPFCGGetTruncationType - get the truncation strategy employed by FCG

449:    Not Collective

451:    Input Parameter:
452: .  ksp - the Krylov space context

454:    Output Parameter:
455: .  truncstrat - the strategy type

457:    Level: intermediate

459: .seealso: KSPFCG, KSPFCGSetTruncationType, KSPFCDTruncationType
460: @*/
461: PetscErrorCode KSPFCGGetTruncationType(KSP ksp,KSPFCDTruncationType *truncstrat)
462: {
463:   KSP_FCG *fcg=(KSP_FCG*)ksp->data;

466:   *truncstrat=fcg->truncstrat;
467:   return 0;
468: }

470: static PetscErrorCode KSPSetFromOptions_FCG(PetscOptionItems *PetscOptionsObject,KSP ksp)
471: {
472:   KSP_FCG        *fcg=(KSP_FCG*)ksp->data;
473:   PetscInt       mmax,nprealloc;
474:   PetscBool      flg;

476:   PetscOptionsHead(PetscOptionsObject,"KSP FCG Options");
477:   PetscOptionsInt("-ksp_fcg_mmax","Maximum number of search directions to store","KSPFCGSetMmax",fcg->mmax,&mmax,&flg);
478:   if (flg) {
479:     KSPFCGSetMmax(ksp,mmax);
480:   }
481:   PetscOptionsInt("-ksp_fcg_nprealloc","Number of directions to preallocate","KSPFCGSetNprealloc",fcg->nprealloc,&nprealloc,&flg);
482:   if (flg) {
483:     KSPFCGSetNprealloc(ksp,nprealloc);
484:   }
485:   PetscOptionsEnum("-ksp_fcg_truncation_type","Truncation approach for directions","KSPFCGSetTruncationType",KSPFCDTruncationTypes,(PetscEnum)fcg->truncstrat,(PetscEnum*)&fcg->truncstrat,NULL);
486:   PetscOptionsTail();
487:   return 0;
488: }

490: /*MC
491:       KSPFCG - Implements the Flexible Conjugate Gradient method (FCG)

493:   Options Database Keys:
494: +   -ksp_fcg_mmax <N>  - maximum number of search directions
495: .   -ksp_fcg_nprealloc <N> - number of directions to preallocate
496: -   -ksp_fcg_truncation_type <standard,notay> - truncation approach for directions

498:     Contributed by Patrick Sanan

500:    Notes:
501:    Supports left preconditioning only.

503:    Level: beginner

505:   References:
506: + * - Notay, Y."Flexible Conjugate Gradients", SIAM J. Sci. Comput. 22:4, 2000
507: - * - Axelsson, O. and Vassilevski, P. S. "A Black Box Generalized Conjugate Gradient Solver with Inner Iterations and Variable step Preconditioning",
508:     SIAM J. Matrix Anal. Appl. 12:4, 1991

510:  .seealso : KSPGCR, KSPFGMRES, KSPCG, KSPFCGSetMmax(), KSPFCGGetMmax(), KSPFCGSetNprealloc(), KSPFCGGetNprealloc(), KSPFCGSetTruncationType(), KSPFCGGetTruncationType()

512: M*/
513: PETSC_EXTERN PetscErrorCode KSPCreate_FCG(KSP ksp)
514: {
515:   KSP_FCG        *fcg;

517:   PetscNewLog(ksp,&fcg);
518: #if !defined(PETSC_USE_COMPLEX)
519:   fcg->type       = KSP_CG_SYMMETRIC;
520: #else
521:   fcg->type       = KSP_CG_HERMITIAN;
522: #endif
523:   fcg->mmax       = KSPFCG_DEFAULT_MMAX;
524:   fcg->nprealloc  = KSPFCG_DEFAULT_NPREALLOC;
525:   fcg->nvecs      = 0;
526:   fcg->vecb       = KSPFCG_DEFAULT_VECB;
527:   fcg->nchunks    = 0;
528:   fcg->truncstrat = KSPFCG_DEFAULT_TRUNCSTRAT;

530:   ksp->data = (void*)fcg;

532:   KSPSetSupportedNorm(ksp,KSP_NORM_PRECONDITIONED,PC_LEFT,2);
533:   KSPSetSupportedNorm(ksp,KSP_NORM_UNPRECONDITIONED,PC_LEFT,1);
534:   KSPSetSupportedNorm(ksp,KSP_NORM_NATURAL,PC_LEFT,1);
535:   KSPSetSupportedNorm(ksp,KSP_NORM_NONE,PC_LEFT,1);

537:   ksp->ops->setup          = KSPSetUp_FCG;
538:   ksp->ops->solve          = KSPSolve_FCG;
539:   ksp->ops->destroy        = KSPDestroy_FCG;
540:   ksp->ops->view           = KSPView_FCG;
541:   ksp->ops->setfromoptions = KSPSetFromOptions_FCG;
542:   ksp->ops->buildsolution  = KSPBuildSolutionDefault;
543:   ksp->ops->buildresidual  = KSPBuildResidualDefault;
544:   return 0;
545: }