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:     KSPCreateVecs(ksp, nnewvecs, &fcg->pPvecs[fcg->nchunks], 0, NULL);
 27:     fcg->nvecs += nnewvecs;
 28:     for (i = 0; i < nnewvecs; ++i) {
 29:       fcg->Cvecs[nvecsprev + i] = fcg->pCvecs[fcg->nchunks][i];
 30:       fcg->Pvecs[nvecsprev + i] = fcg->pPvecs[fcg->nchunks][i];
 31:     }
 32:     fcg->chunksizes[fcg->nchunks] = nnewvecs;
 33:     ++fcg->nchunks;
 34:   }
 35:   return 0;
 36: }

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


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

 48:   /* Allocated space for pointers to additional work vectors
 49:    note that mmax is the number of previous directions, so we add 1 for the current direction,
 50:    and an extra 1 for the prealloc (which might be empty) */
 51:   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);

 53:   /* If the requested number of preallocated vectors is greater than mmax reduce nprealloc */
 54:   if (fcg->nprealloc > fcg->mmax + 1) PetscInfo(NULL, "Requested nprealloc=%" PetscInt_FMT " is greater than m_max+1=%" PetscInt_FMT ". Resetting nprealloc = m_max+1.\n", fcg->nprealloc, fcg->mmax + 1);

 56:   /* Preallocate additional work vectors */
 57:   KSPAllocateVectors_FCG(ksp, fcg->nprealloc, fcg->nprealloc);
 58:   /*
 59:   If user requested computations of eigenvalues then allocate work
 60:   work space needed
 61:   */
 62:   if (ksp->calc_sings) {
 63:     /* get space to store tridiagonal matrix for Lanczos */
 64:     PetscMalloc4(maxit, &fcg->e, maxit, &fcg->d, maxit, &fcg->ee, maxit, &fcg->dd);

 66:     ksp->ops->computeextremesingularvalues = KSPComputeExtremeSingularValues_CG;
 67:     ksp->ops->computeeigenvalues           = KSPComputeEigenvalues_CG;
 68:   }
 69:   return 0;
 70: }

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


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

 88:   X = ksp->vec_sol;
 89:   B = ksp->vec_rhs;
 90:   R = ksp->work[0];
 91:   Z = ksp->work[1];

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

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

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

144:   i = 0;
145:   do {
146:     ksp->its = i + 1;

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

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

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

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

171:     {
172:       PetscInt l, ndots;

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

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

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

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

229:     /* Check for convergence */
230:     ksp->rnorm = dp;
231:     KSPLogResidualHistory(ksp, dp);
232:     KSPMonitor(ksp, i + 1, dp);
233:     (*ksp->converged)(ksp, i + 1, dp, &ksp->reason, ksp->cnvP);
234:     if (ksp->reason) break;

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

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

242:     if (eigs) {
243:       if (i > 0) {
245:         e[i] = PetscSqrtReal(PetscAbsScalar(beta / betaold)) / alphaold;
246:         d[i] = PetscSqrtReal(PetscAbsScalar(beta / betaold)) * e[i] + 1.0 / alpha;
247:       } else {
248:         d[i] = PetscSqrtReal(PetscAbsScalar(beta)) * e[i] + 1.0 / alpha;
249:       }
250:       fcg->ned = ksp->its - 1;
251:     }
252:     ++i;
253:   } while (i < ksp->max_it);
254:   if (i >= ksp->max_it) ksp->reason = KSP_DIVERGED_ITS;
255:   return 0;
256: }

258: static PetscErrorCode KSPDestroy_FCG(KSP ksp)
259: {
260:   PetscInt i;
261:   KSP_FCG *fcg = (KSP_FCG *)ksp->data;


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

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

281: static PetscErrorCode KSPView_FCG(KSP ksp, PetscViewer viewer)
282: {
283:   KSP_FCG    *fcg = (KSP_FCG *)ksp->data;
284:   PetscBool   iascii, isstring;
285:   const char *truncstr;

287:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
288:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSTRING, &isstring);

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

294:   if (iascii) {
295:     PetscViewerASCIIPrintf(viewer, "  m_max=%" PetscInt_FMT "\n", fcg->mmax);
296:     PetscViewerASCIIPrintf(viewer, "  preallocated %" PetscInt_FMT " directions\n", PetscMin(fcg->nprealloc, fcg->mmax + 1));
297:     PetscViewerASCIIPrintf(viewer, "  %s\n", truncstr);
298:   } else if (isstring) {
299:     PetscViewerStringSPrintf(viewer, "m_max %" PetscInt_FMT " nprealloc %" PetscInt_FMT " %s", fcg->mmax, fcg->nprealloc, truncstr);
300:   }
301:   return 0;
302: }

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

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

311:   Logically Collective on ksp

313:   Input Parameters:
314: +  ksp - the Krylov space context
315: -  mmax - the maximum number of previous directions to orthogonalize againt

317:   Level: intermediate

319: .seealso: `KSPFCG`, `KSPFCGGetTruncationType()`, `KSPFCGGetNprealloc()`
320: @*/
321: PetscErrorCode KSPFCGSetMmax(KSP ksp, PetscInt mmax)
322: {
323:   KSP_FCG *fcg = (KSP_FCG *)ksp->data;

327:   fcg->mmax = mmax;
328:   return 0;
329: }

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

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

336:    Not Collective

338:    Input Parameter:
339: .  ksp - the Krylov space context

341:    Output Parameter:
342: .  mmax - the maximum number of previous directions allowed for orthogonalization

344:    Level: intermediate

346: .seealso: `KSPFCG`, `KSPFCGGetTruncationType()`, `KSPFCGGetNprealloc()`, `KSPFCGSetMmax()`
347: @*/

349: PetscErrorCode KSPFCGGetMmax(KSP ksp, PetscInt *mmax)
350: {
351:   KSP_FCG *fcg = (KSP_FCG *)ksp->data;

354:   *mmax = fcg->mmax;
355:   return 0;
356: }

358: /*@
359:   KSPFCGSetNprealloc - set the number of directions to preallocate with FCG

361:   Logically Collective on ksp

363:   Input Parameters:
364: +  ksp - the Krylov space context
365: -  nprealloc - the number of vectors to preallocate

367:   Level: advanced

369:   Options Database:
370: . -ksp_fcg_nprealloc <N> - number of directions to preallocate

372: .seealso: `KSPFCG`, `KSPFCGGetTruncationType()`, `KSPFCGGetNprealloc()`
373: @*/
374: PetscErrorCode KSPFCGSetNprealloc(KSP ksp, PetscInt nprealloc)
375: {
376:   KSP_FCG *fcg = (KSP_FCG *)ksp->data;

381:   fcg->nprealloc = nprealloc;
382:   return 0;
383: }

385: /*@
386:   KSPFCGGetNprealloc - get the number of directions preallocate by FCG

388:    Not Collective

390:    Input Parameter:
391: .  ksp - the Krylov space context

393:    Output Parameter:
394: .  nprealloc - the number of directions preallocated

396:    Level: advanced

398: .seealso: `KSPFCG`, `KSPFCGGetTruncationType()`, `KSPFCGSetNprealloc()`
399: @*/
400: PetscErrorCode KSPFCGGetNprealloc(KSP ksp, PetscInt *nprealloc)
401: {
402:   KSP_FCG *fcg = (KSP_FCG *)ksp->data;

405:   *nprealloc = fcg->nprealloc;
406:   return 0;
407: }

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

412:   Logically Collective on ksp

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

417:   Input Parameters:
418: +  ksp - the Krylov space context
419: -  truncstrat - the choice of strategy

421:   Level: intermediate

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

426:   .seealso: `KSPFCDTruncationType`, `KSPFCGGetTruncationType`
427: @*/
428: PetscErrorCode KSPFCGSetTruncationType(KSP ksp, KSPFCDTruncationType truncstrat)
429: {
430:   KSP_FCG *fcg = (KSP_FCG *)ksp->data;

434:   fcg->truncstrat = truncstrat;
435:   return 0;
436: }

438: /*@
439:   KSPFCGGetTruncationType - get the truncation strategy employed by FCG

441:    Not Collective

443:    Input Parameter:
444: .  ksp - the Krylov space context

446:    Output Parameter:
447: .  truncstrat - the strategy type

449:    Level: intermediate

451: .seealso: `KSPFCG`, `KSPFCGSetTruncationType`, `KSPFCDTruncationType`
452: @*/
453: PetscErrorCode KSPFCGGetTruncationType(KSP ksp, KSPFCDTruncationType *truncstrat)
454: {
455:   KSP_FCG *fcg = (KSP_FCG *)ksp->data;

458:   *truncstrat = fcg->truncstrat;
459:   return 0;
460: }

462: static PetscErrorCode KSPSetFromOptions_FCG(KSP ksp, PetscOptionItems *PetscOptionsObject)
463: {
464:   KSP_FCG  *fcg = (KSP_FCG *)ksp->data;
465:   PetscInt  mmax, nprealloc;
466:   PetscBool flg;

468:   PetscOptionsHeadBegin(PetscOptionsObject, "KSP FCG Options");
469:   PetscOptionsInt("-ksp_fcg_mmax", "Maximum number of search directions to store", "KSPFCGSetMmax", fcg->mmax, &mmax, &flg);
470:   if (flg) KSPFCGSetMmax(ksp, mmax);
471:   PetscOptionsInt("-ksp_fcg_nprealloc", "Number of directions to preallocate", "KSPFCGSetNprealloc", fcg->nprealloc, &nprealloc, &flg);
472:   if (flg) KSPFCGSetNprealloc(ksp, nprealloc);
473:   PetscOptionsEnum("-ksp_fcg_truncation_type", "Truncation approach for directions", "KSPFCGSetTruncationType", KSPFCDTruncationTypes, (PetscEnum)fcg->truncstrat, (PetscEnum *)&fcg->truncstrat, NULL);
474:   PetscOptionsHeadEnd();
475:   return 0;
476: }

478: /*MC
479:       KSPFCG - Implements the Flexible Conjugate Gradient method (FCG)

481:   Options Database Keys:
482: +   -ksp_fcg_mmax <N>  - maximum number of search directions
483: .   -ksp_fcg_nprealloc <N> - number of directions to preallocate
484: -   -ksp_fcg_truncation_type <standard,notay> - truncation approach for directions

486:     Contributed by Patrick Sanan

488:    Notes:
489:    Supports left preconditioning only.

491:    Level: beginner

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

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

500: M*/
501: PETSC_EXTERN PetscErrorCode KSPCreate_FCG(KSP ksp)
502: {
503:   KSP_FCG *fcg;

505:   PetscNew(&fcg);
506: #if !defined(PETSC_USE_COMPLEX)
507:   fcg->type = KSP_CG_SYMMETRIC;
508: #else
509:   fcg->type = KSP_CG_HERMITIAN;
510: #endif
511:   fcg->mmax       = KSPFCG_DEFAULT_MMAX;
512:   fcg->nprealloc  = KSPFCG_DEFAULT_NPREALLOC;
513:   fcg->nvecs      = 0;
514:   fcg->vecb       = KSPFCG_DEFAULT_VECB;
515:   fcg->nchunks    = 0;
516:   fcg->truncstrat = KSPFCG_DEFAULT_TRUNCSTRAT;

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

520:   KSPSetSupportedNorm(ksp, KSP_NORM_PRECONDITIONED, PC_LEFT, 2);
521:   KSPSetSupportedNorm(ksp, KSP_NORM_UNPRECONDITIONED, PC_LEFT, 1);
522:   KSPSetSupportedNorm(ksp, KSP_NORM_NATURAL, PC_LEFT, 1);
523:   KSPSetSupportedNorm(ksp, KSP_NORM_NONE, PC_LEFT, 1);

525:   ksp->ops->setup          = KSPSetUp_FCG;
526:   ksp->ops->solve          = KSPSolve_FCG;
527:   ksp->ops->destroy        = KSPDestroy_FCG;
528:   ksp->ops->view           = KSPView_FCG;
529:   ksp->ops->setfromoptions = KSPSetFromOptions_FCG;
530:   ksp->ops->buildsolution  = KSPBuildSolutionDefault;
531:   ksp->ops->buildresidual  = KSPBuildResidualDefault;
532:   return 0;
533: }