Actual source code: rii.c

slepc-3.12.0 2019-09-30
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    SLEPc nonlinear eigensolver: "rii"

 13:    Method: Residual inverse iteration

 15:    Algorithm:

 17:        Simple residual inverse iteration with varying shift.

 19:    References:

 21:        [1] A. Neumaier, "Residual inverse iteration for the nonlinear
 22:            eigenvalue problem", SIAM J. Numer. Anal. 22(5):914-923, 1985.
 23: */

 25: #include <slepc/private/nepimpl.h>         /*I "slepcnep.h" I*/
 26: #include <../src/nep/impls/nepdefl.h>

 28: typedef struct {
 29:   PetscInt  max_inner_it;     /* maximum number of Newton iterations */
 30:   PetscInt  lag;              /* interval to rebuild preconditioner */
 31:   PetscBool cctol;            /* constant correction tolerance */
 32:   KSP       ksp;              /* linear solver object */
 33: } NEP_RII;

 35: PetscErrorCode NEPSetUp_RII(NEP nep)
 36: {
 38:   PetscBool      istrivial;

 41:   if (nep->ncv) { PetscInfo(nep,"Setting ncv = nev, ignoring user-provided value\n"); }
 42:   nep->ncv = nep->nev;
 43:   if (nep->mpd) { PetscInfo(nep,"Setting mpd = nev, ignoring user-provided value\n"); }
 44:   nep->mpd = nep->nev;
 45:   if (!nep->max_it) nep->max_it = PetscMax(5000,2*nep->n/nep->ncv);
 46:   if (nep->which && nep->which!=NEP_TARGET_MAGNITUDE) SETERRQ(PetscObjectComm((PetscObject)nep),1,"Wrong value of which");

 48:   RGIsTrivial(nep->rg,&istrivial);
 49:   if (!istrivial) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"This solver does not support region filtering");

 51:   NEPAllocateSolution(nep,0);
 52:   NEPSetWorkVecs(nep,2);
 53:   return(0);
 54: }

 56: PetscErrorCode NEPSolve_RII(NEP nep)
 57: {
 58:   PetscErrorCode     ierr;
 59:   NEP_RII            *ctx = (NEP_RII*)nep->data;
 60:   Mat                T,Tp,H;
 61:   Vec                uu,u,r,delta;
 62:   PetscScalar        lambda,lambda2,sigma,a1,a2,corr,*Hp,*Ap;
 63:   PetscReal          nrm,resnorm=1.0,ktol=0.1,perr;
 64:   PetscBool          skip=PETSC_FALSE;
 65:   PetscInt           inner_its,its=0,ldh,ldds,i,j;
 66:   NEP_EXT_OP         extop=NULL;
 67:   KSPConvergedReason kspreason;

 70:   /* get initial approximation of eigenvalue and eigenvector */
 71:   NEPGetDefaultShift(nep,&sigma);
 72:   lambda = sigma;
 73:   if (!nep->nini) {
 74:     BVSetRandomColumn(nep->V,0);
 75:     BVNormColumn(nep->V,0,NORM_2,&nrm);
 76:     BVScaleColumn(nep->V,0,1.0/nrm);
 77:   }
 78:   if (!ctx->ksp) { NEPRIIGetKSP(nep,&ctx->ksp); }
 79:   NEPDeflationInitialize(nep,nep->V,ctx->ksp,PETSC_FALSE,nep->nev,&extop);
 80:   NEPDeflationCreateVec(extop,&u);
 81:   VecDuplicate(u,&r);
 82:   VecDuplicate(u,&delta);
 83:   BVGetColumn(nep->V,0,&uu);
 84:   NEPDeflationCopyToExtendedVec(extop,uu,NULL,u,PETSC_FALSE);
 85:   BVRestoreColumn(nep->V,0,&uu);

 87:   /* prepare linear solver */
 88:   NEPDeflationSolveSetUp(extop,sigma);

 90:   VecCopy(u,r);
 91:   NEPDeflationFunctionSolve(extop,r,u);
 92:   VecNorm(u,NORM_2,&nrm);
 93:   VecScale(u,1.0/nrm);

 95:   /* Restart loop */
 96:   while (nep->reason == NEP_CONVERGED_ITERATING) {
 97:     its++;

 99:     /* Use Newton's method to compute nonlinear Rayleigh functional. Current eigenvalue
100:        estimate as starting value. */
101:     inner_its=0;
102:     lambda2 = lambda;
103:     do {
104:       NEPDeflationComputeFunction(extop,lambda,&T);
105:       MatMult(T,u,r);
106:       VecDot(r,u,&a1);
107:       if (inner_its && PetscAbsScalar(a1)/PetscAbsScalar(lambda)<=PETSC_SQRT_MACHINE_EPSILON) break;
108:       NEPDeflationComputeJacobian(extop,lambda,&Tp);
109:       MatMult(Tp,u,r);
110:       VecDot(r,u,&a2);
111:       corr = a1/a2;
112:       lambda = lambda - corr;
113:       inner_its++;
114:     } while (PetscAbsScalar(corr)/PetscAbsScalar(lambda)>PETSC_SQRT_MACHINE_EPSILON && inner_its<ctx->max_inner_it);

116:     /* form residual,  r = T(lambda)*u */
117:     NEPDeflationComputeFunction(extop,lambda,&T);
118:     MatMult(T,u,r);

120:     /* convergence test */
121:     perr = nep->errest[nep->nconv];
122:     VecNorm(r,NORM_2,&resnorm);
123:     (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
124:     nep->eigr[nep->nconv] = lambda;
125:     if (nep->errest[nep->nconv]<=nep->tol) {
126:       nep->nconv = nep->nconv + 1;
127:       NEPDeflationLocking(extop,u,lambda);
128:       nep->its += its;
129:       skip = PETSC_TRUE;
130:     }
131:     (*nep->stopping)(nep,nep->its+its,nep->max_it,nep->nconv,nep->nev,&nep->reason,nep->stoppingctx);
132:     if (!skip || nep->reason>0) {
133:       NEPMonitor(nep,nep->its+its,nep->nconv,nep->eigr,nep->eigi,nep->errest,(nep->reason>0)?nep->nconv:nep->nconv+1);
134:     }

136:     if (nep->reason == NEP_CONVERGED_ITERATING) {
137:       if (!skip) {
138:         /* update preconditioner and set adaptive tolerance */
139:         if (ctx->lag && !(its%ctx->lag) && its>=2*ctx->lag && perr && nep->errest[nep->nconv]>.5*perr) {
140:           NEPDeflationSolveSetUp(extop,lambda2);
141:         }
142:         if (!ctx->cctol) {
143:           ktol = PetscMax(ktol/2.0,PETSC_MACHINE_EPSILON*10.0);
144:           KSPSetTolerances(ctx->ksp,ktol,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
145:         }

147:         /* eigenvector correction: delta = T(sigma)\r */
148:         NEPDeflationFunctionSolve(extop,r,delta);
149:         KSPGetConvergedReason(ctx->ksp,&kspreason);
150:         if (kspreason<0) {
151:           PetscInfo1(nep,"iter=%D, linear solve failed, stopping solve\n",nep->its);
152:           nep->reason = NEP_DIVERGED_LINEAR_SOLVE;
153:           break;
154:         }

156:         /* update eigenvector: u = u - delta */
157:         VecAXPY(u,-1.0,delta);

159:         /* normalize eigenvector */
160:         VecNormalize(u,NULL);
161:       } else {
162:         its = -1;
163:         NEPDeflationSetRandomVec(extop,u);
164:         NEPDeflationSolveSetUp(extop,sigma);
165:         VecCopy(u,r);
166:         NEPDeflationFunctionSolve(extop,r,u);
167:         VecNorm(u,NORM_2,&nrm);
168:         VecScale(u,1.0/nrm);
169:         lambda = sigma;
170:         skip = PETSC_FALSE;
171:       }
172:     }
173:   }
174:   NEPDeflationGetInvariantPair(extop,NULL,&H);
175:   MatGetSize(H,NULL,&ldh);
176:   DSSetType(nep->ds,DSNHEP);
177:   DSAllocate(nep->ds,PetscMax(nep->nconv,1));
178:   DSGetLeadingDimension(nep->ds,&ldds);
179:   MatDenseGetArray(H,&Hp);
180:   DSGetArray(nep->ds,DS_MAT_A,&Ap);
181:   for (j=0;j<nep->nconv;j++)
182:     for (i=0;i<nep->nconv;i++) Ap[j*ldds+i] = Hp[j*ldh+i];
183:   DSRestoreArray(nep->ds,DS_MAT_A,&Ap);
184:   MatDenseRestoreArray(H,&Hp);
185:   MatDestroy(&H);
186:   DSSetDimensions(nep->ds,nep->nconv,0,0,nep->nconv);
187:   DSSolve(nep->ds,nep->eigr,nep->eigi);
188:   NEPDeflationReset(extop);
189:   VecDestroy(&u);
190:   VecDestroy(&r);
191:   VecDestroy(&delta);
192:   return(0);
193: }

195: PetscErrorCode NEPSetFromOptions_RII(PetscOptionItems *PetscOptionsObject,NEP nep)
196: {
198:   NEP_RII        *ctx = (NEP_RII*)nep->data;
199:   PetscBool      flg;
200:   PetscInt       i;

203:   PetscOptionsHead(PetscOptionsObject,"NEP RII Options");

205:     i = 0;
206:     PetscOptionsInt("-nep_rii_max_it","Maximum number of Newton iterations for updating Rayleigh functional","NEPRIISetMaximumIterations",ctx->max_inner_it,&i,&flg);
207:     if (flg) { NEPRIISetMaximumIterations(nep,i); }

209:     PetscOptionsBool("-nep_rii_const_correction_tol","Constant correction tolerance for the linear solver","NEPRIISetConstCorrectionTol",ctx->cctol,&ctx->cctol,NULL);

211:     i = 0;
212:     PetscOptionsInt("-nep_rii_lag_preconditioner","Interval to rebuild preconditioner","NEPRIISetLagPreconditioner",ctx->lag,&i,&flg);
213:     if (flg) { NEPRIISetLagPreconditioner(nep,i); }

215:   PetscOptionsTail();

217:   if (!ctx->ksp) { NEPRIIGetKSP(nep,&ctx->ksp); }
218:   KSPSetFromOptions(ctx->ksp);
219:   return(0);
220: }

222: static PetscErrorCode NEPRIISetMaximumIterations_RII(NEP nep,PetscInt its)
223: {
224:   NEP_RII *ctx = (NEP_RII*)nep->data;

227:   if (its==PETSC_DEFAULT) ctx->max_inner_it = 10;
228:   else {
229:     if (its<=0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of iterations must be >0");
230:     ctx->max_inner_it = its;
231:   }
232:   return(0);
233: }

235: /*@
236:    NEPRIISetMaximumIterations - Sets the maximum number of inner iterations to be
237:    used in the RII solver. These are the Newton iterations related to the computation
238:    of the nonlinear Rayleigh functional.

240:    Logically Collective on nep

242:    Input Parameters:
243: +  nep - nonlinear eigenvalue solver
244: -  its - maximum inner iterations

246:    Level: advanced

248: .seealso: NEPRIIGetMaximumIterations()
249: @*/
250: PetscErrorCode NEPRIISetMaximumIterations(NEP nep,PetscInt its)
251: {

257:   PetscTryMethod(nep,"NEPRIISetMaximumIterations_C",(NEP,PetscInt),(nep,its));
258:   return(0);
259: }

261: static PetscErrorCode NEPRIIGetMaximumIterations_RII(NEP nep,PetscInt *its)
262: {
263:   NEP_RII *ctx = (NEP_RII*)nep->data;

266:   *its = ctx->max_inner_it;
267:   return(0);
268: }

270: /*@
271:    NEPRIIGetMaximumIterations - Gets the maximum number of inner iterations of RII.

273:    Not Collective

275:    Input Parameter:
276: .  nep - nonlinear eigenvalue solver

278:    Output Parameter:
279: .  its - maximum inner iterations

281:    Level: advanced

283: .seealso: NEPRIISetMaximumIterations()
284: @*/
285: PetscErrorCode NEPRIIGetMaximumIterations(NEP nep,PetscInt *its)
286: {

292:   PetscUseMethod(nep,"NEPRIIGetMaximumIterations_C",(NEP,PetscInt*),(nep,its));
293:   return(0);
294: }

296: static PetscErrorCode NEPRIISetLagPreconditioner_RII(NEP nep,PetscInt lag)
297: {
298:   NEP_RII *ctx = (NEP_RII*)nep->data;

301:   if (lag<0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag must be non-negative");
302:   ctx->lag = lag;
303:   return(0);
304: }

306: /*@
307:    NEPRIISetLagPreconditioner - Determines when the preconditioner is rebuilt in the
308:    nonlinear solve.

310:    Logically Collective on nep

312:    Input Parameters:
313: +  nep - nonlinear eigenvalue solver
314: -  lag - 0 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is
315:           computed within the nonlinear iteration, 2 means every second time
316:           the Jacobian is built, etc.

318:    Options Database Keys:
319: .  -nep_rii_lag_preconditioner <lag>

321:    Notes:
322:    The default is 1.
323:    The preconditioner is ALWAYS built in the first iteration of a nonlinear solve.

325:    Level: intermediate

327: .seealso: NEPRIIGetLagPreconditioner()
328: @*/
329: PetscErrorCode NEPRIISetLagPreconditioner(NEP nep,PetscInt lag)
330: {

336:   PetscTryMethod(nep,"NEPRIISetLagPreconditioner_C",(NEP,PetscInt),(nep,lag));
337:   return(0);
338: }

340: static PetscErrorCode NEPRIIGetLagPreconditioner_RII(NEP nep,PetscInt *lag)
341: {
342:   NEP_RII *ctx = (NEP_RII*)nep->data;

345:   *lag = ctx->lag;
346:   return(0);
347: }

349: /*@
350:    NEPRIIGetLagPreconditioner - Indicates how often the preconditioner is rebuilt.

352:    Not Collective

354:    Input Parameter:
355: .  nep - nonlinear eigenvalue solver

357:    Output Parameter:
358: .  lag - the lag parameter

360:    Level: intermediate

362: .seealso: NEPRIISetLagPreconditioner()
363: @*/
364: PetscErrorCode NEPRIIGetLagPreconditioner(NEP nep,PetscInt *lag)
365: {

371:   PetscUseMethod(nep,"NEPRIIGetLagPreconditioner_C",(NEP,PetscInt*),(nep,lag));
372:   return(0);
373: }

375: static PetscErrorCode NEPRIISetConstCorrectionTol_RII(NEP nep,PetscBool cct)
376: {
377:   NEP_RII *ctx = (NEP_RII*)nep->data;

380:   ctx->cctol = cct;
381:   return(0);
382: }

384: /*@
385:    NEPRIISetConstCorrectionTol - Sets a flag to keep the tolerance used
386:    in the linear solver constant.

388:    Logically Collective on nep

390:    Input Parameters:
391: +  nep - nonlinear eigenvalue solver
392: -  cct - a boolean value

394:    Options Database Keys:
395: .  -nep_rii_const_correction_tol <bool> - set the boolean flag

397:    Notes:
398:    By default, an exponentially decreasing tolerance is set in the KSP used
399:    within the nonlinear iteration, so that each Newton iteration requests
400:    better accuracy than the previous one. The constant correction tolerance
401:    flag stops this behaviour.

403:    Level: intermediate

405: .seealso: NEPRIIGetConstCorrectionTol()
406: @*/
407: PetscErrorCode NEPRIISetConstCorrectionTol(NEP nep,PetscBool cct)
408: {

414:   PetscTryMethod(nep,"NEPRIISetConstCorrectionTol_C",(NEP,PetscBool),(nep,cct));
415:   return(0);
416: }

418: static PetscErrorCode NEPRIIGetConstCorrectionTol_RII(NEP nep,PetscBool *cct)
419: {
420:   NEP_RII *ctx = (NEP_RII*)nep->data;

423:   *cct = ctx->cctol;
424:   return(0);
425: }

427: /*@
428:    NEPRIIGetConstCorrectionTol - Returns the constant tolerance flag.

430:    Not Collective

432:    Input Parameter:
433: .  nep - nonlinear eigenvalue solver

435:    Output Parameter:
436: .  cct - the value of the constant tolerance flag

438:    Level: intermediate

440: .seealso: NEPRIISetConstCorrectionTol()
441: @*/
442: PetscErrorCode NEPRIIGetConstCorrectionTol(NEP nep,PetscBool *cct)
443: {

449:   PetscUseMethod(nep,"NEPRIIGetConstCorrectionTol_C",(NEP,PetscBool*),(nep,cct));
450:   return(0);
451: }

453: static PetscErrorCode NEPRIISetKSP_RII(NEP nep,KSP ksp)
454: {
456:   NEP_RII        *ctx = (NEP_RII*)nep->data;

459:   PetscObjectReference((PetscObject)ksp);
460:   KSPDestroy(&ctx->ksp);
461:   ctx->ksp = ksp;
462:   PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->ksp);
463:   nep->state = NEP_STATE_INITIAL;
464:   return(0);
465: }

467: /*@
468:    NEPRIISetKSP - Associate a linear solver object (KSP) to the nonlinear
469:    eigenvalue solver.

471:    Collective on nep

473:    Input Parameters:
474: +  nep - eigenvalue solver
475: -  ksp - the linear solver object

477:    Level: advanced

479: .seealso: NEPRIIGetKSP()
480: @*/
481: PetscErrorCode NEPRIISetKSP(NEP nep,KSP ksp)
482: {

489:   PetscTryMethod(nep,"NEPRIISetKSP_C",(NEP,KSP),(nep,ksp));
490:   return(0);
491: }

493: static PetscErrorCode NEPRIIGetKSP_RII(NEP nep,KSP *ksp)
494: {
496:   NEP_RII        *ctx = (NEP_RII*)nep->data;

499:   if (!ctx->ksp) {
500:     KSPCreate(PetscObjectComm((PetscObject)nep),&ctx->ksp);
501:     PetscObjectIncrementTabLevel((PetscObject)ctx->ksp,(PetscObject)nep,1);
502:     KSPSetOptionsPrefix(ctx->ksp,((PetscObject)nep)->prefix);
503:     KSPAppendOptionsPrefix(ctx->ksp,"nep_rii_");
504:     PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->ksp);
505:     PetscObjectSetOptions((PetscObject)ctx->ksp,((PetscObject)nep)->options);
506:     KSPSetErrorIfNotConverged(ctx->ksp,PETSC_TRUE);
507:     KSPSetTolerances(ctx->ksp,SLEPC_DEFAULT_TOL,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
508:   }
509:   *ksp = ctx->ksp;
510:   return(0);
511: }

513: /*@
514:    NEPRIIGetKSP - Retrieve the linear solver object (KSP) associated with
515:    the nonlinear eigenvalue solver.

517:    Not Collective

519:    Input Parameter:
520: .  nep - nonlinear eigenvalue solver

522:    Output Parameter:
523: .  ksp - the linear solver object

525:    Level: advanced

527: .seealso: NEPRIISetKSP()
528: @*/
529: PetscErrorCode NEPRIIGetKSP(NEP nep,KSP *ksp)
530: {

536:   PetscUseMethod(nep,"NEPRIIGetKSP_C",(NEP,KSP*),(nep,ksp));
537:   return(0);
538: }

540: PetscErrorCode NEPView_RII(NEP nep,PetscViewer viewer)
541: {
543:   NEP_RII        *ctx = (NEP_RII*)nep->data;
544:   PetscBool      isascii;

547:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
548:   if (isascii) {
549:     PetscViewerASCIIPrintf(viewer,"  maximum number of inner iterations: %D\n",ctx->max_inner_it);
550:     if (ctx->cctol) {
551:       PetscViewerASCIIPrintf(viewer,"  using a constant tolerance for the linear solver\n");
552:     }
553:     if (ctx->lag) {
554:       PetscViewerASCIIPrintf(viewer,"  updating the preconditioner every %D iterations\n",ctx->lag);
555:     }
556:     if (!ctx->ksp) { NEPRIIGetKSP(nep,&ctx->ksp); }
557:     PetscViewerASCIIPushTab(viewer);
558:     KSPView(ctx->ksp,viewer);
559:     PetscViewerASCIIPopTab(viewer);
560:   }
561:   return(0);
562: }

564: PetscErrorCode NEPReset_RII(NEP nep)
565: {
567:   NEP_RII        *ctx = (NEP_RII*)nep->data;

570:   KSPReset(ctx->ksp);
571:   return(0);
572: }

574: PetscErrorCode NEPDestroy_RII(NEP nep)
575: {
577:   NEP_RII        *ctx = (NEP_RII*)nep->data;

580:   KSPDestroy(&ctx->ksp);
581:   PetscFree(nep->data);
582:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetMaximumIterations_C",NULL);
583:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetMaximumIterations_C",NULL);
584:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetLagPreconditioner_C",NULL);
585:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetLagPreconditioner_C",NULL);
586:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetConstCorrectionTol_C",NULL);
587:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetConstCorrectionTol_C",NULL);
588:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetKSP_C",NULL);
589:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetKSP_C",NULL);
590:   return(0);
591: }

593: SLEPC_EXTERN PetscErrorCode NEPCreate_RII(NEP nep)
594: {
596:   NEP_RII        *ctx;

599:   PetscNewLog(nep,&ctx);
600:   nep->data = (void*)ctx;
601:   ctx->max_inner_it = 10;
602:   ctx->lag          = 1;
603:   ctx->cctol        = PETSC_FALSE;

605:   nep->useds = PETSC_TRUE;

607:   nep->ops->solve          = NEPSolve_RII;
608:   nep->ops->setup          = NEPSetUp_RII;
609:   nep->ops->setfromoptions = NEPSetFromOptions_RII;
610:   nep->ops->reset          = NEPReset_RII;
611:   nep->ops->destroy        = NEPDestroy_RII;
612:   nep->ops->view           = NEPView_RII;
613:   nep->ops->computevectors = NEPComputeVectors_Schur;

615:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetMaximumIterations_C",NEPRIISetMaximumIterations_RII);
616:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetMaximumIterations_C",NEPRIIGetMaximumIterations_RII);
617:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetLagPreconditioner_C",NEPRIISetLagPreconditioner_RII);
618:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetLagPreconditioner_C",NEPRIIGetLagPreconditioner_RII);
619:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetConstCorrectionTol_C",NEPRIISetConstCorrectionTol_RII);
620:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetConstCorrectionTol_C",NEPRIIGetConstCorrectionTol_RII);
621:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetKSP_C",NEPRIISetKSP_RII);
622:   PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetKSP_C",NEPRIIGetKSP_RII);
623:   return(0);
624: }