Actual source code: qarnoldi.c

slepc-3.7.2 2016-07-19
Report Typos and Errors
  1: /*

  3:    SLEPc quadratic eigensolver: "qarnoldi"

  5:    Method: Q-Arnoldi

  7:    Algorithm:

  9:        Quadratic Arnoldi with Krylov-Schur type restart.

 11:    References:

 13:        [1] K. Meerbergen, "The Quadratic Arnoldi method for the solution
 14:            of the quadratic eigenvalue problem", SIAM J. Matrix Anal.
 15:            Appl. 30(4):1462-1482, 2008.

 17:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 18:    SLEPc - Scalable Library for Eigenvalue Problem Computations
 19:    Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain

 21:    This file is part of SLEPc.

 23:    SLEPc is free software: you can redistribute it and/or modify it under  the
 24:    terms of version 3 of the GNU Lesser General Public License as published by
 25:    the Free Software Foundation.

 27:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY
 28:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS
 29:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for
 30:    more details.

 32:    You  should have received a copy of the GNU Lesser General  Public  License
 33:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 34:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 35: */

 37: #include <slepc/private/pepimpl.h>    /*I "slepcpep.h" I*/
 38: #include <petscblaslapack.h>

 40: typedef struct {
 41:   PetscReal keep;         /* restart parameter */
 42:   PetscBool lock;         /* locking/non-locking variant */
 43: } PEP_QARNOLDI;

 47: PetscErrorCode PEPSetUp_QArnoldi(PEP pep)
 48: {
 50:   PEP_QARNOLDI   *ctx = (PEP_QARNOLDI*)pep->data;
 51:   PetscBool      shift,sinv,flg;

 54:   pep->lineariz = PETSC_TRUE;
 55:   PEPSetDimensions_Default(pep,pep->nev,&pep->ncv,&pep->mpd);
 56:   if (!ctx->lock && pep->mpd<pep->ncv) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Should not use mpd parameter in non-locking variant");
 57:   if (!pep->max_it) pep->max_it = PetscMax(100,4*pep->n/pep->ncv);
 58:   /* Set STSHIFT as the default ST */
 59:   if (!((PetscObject)pep->st)->type_name) {
 60:     STSetType(pep->st,STSHIFT);
 61:   }
 62:   PetscObjectTypeCompare((PetscObject)pep->st,STSHIFT,&shift);
 63:   PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv);
 64:   if (!shift && !sinv) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Only STSHIFT and STSINVERT spectral transformations can be used");
 65:   if (!pep->which) {
 66:     if (sinv) pep->which = PEP_TARGET_MAGNITUDE;
 67:     else pep->which = PEP_LARGEST_MAGNITUDE;
 68:   }

 70:   if (pep->nmat!=3) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Solver only available for quadratic problems");
 71:   if (pep->basis!=PEP_BASIS_MONOMIAL) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Solver not implemented for non-monomial bases");
 72:   STGetTransform(pep->st,&flg);
 73:   if (!flg) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Solver requires the ST transformation flag set, see STSetTransform()");

 75:   /* set default extraction */
 76:   if (!pep->extract) {
 77:     pep->extract = PEP_EXTRACT_NONE;
 78:   }
 79:   if (pep->extract!=PEP_EXTRACT_NONE) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Solver does not support requested extraction");
 80:  
 81:   if (!ctx->keep) ctx->keep = 0.5;

 83:   PEPAllocateSolution(pep,0);
 84:   PEPSetWorkVecs(pep,4);

 86:   DSSetType(pep->ds,DSNHEP);
 87:   DSSetExtraRow(pep->ds,PETSC_TRUE);
 88:   DSAllocate(pep->ds,pep->ncv+1);

 90:   /* process starting vector */
 91:   if (pep->nini>-2) {
 92:     BVSetRandomColumn(pep->V,0);
 93:     BVSetRandomColumn(pep->V,1);
 94:   } else {
 95:     BVInsertVec(pep->V,0,pep->IS[0]);
 96:     BVInsertVec(pep->V,1,pep->IS[1]);
 97:   }
 98:   if (pep->nini<0) {
 99:     SlepcBasisDestroy_Private(&pep->nini,&pep->IS);
100:   }
101:   return(0);
102: }

106: PetscErrorCode PEPExtractVectors_QArnoldi(PEP pep)
107: {
109:   PetscInt       i,k=pep->nconv,ldds;
110:   PetscScalar    *X,*pX0;
111:   Mat            X0;

114:   if (pep->nconv==0) return(0);
115:   DSGetLeadingDimension(pep->ds,&ldds);
116:   DSVectors(pep->ds,DS_MAT_X,NULL,NULL);
117:   DSGetArray(pep->ds,DS_MAT_X,&X);

119:   /* update vectors V = V*X */ 
120:   MatCreateSeqDense(PETSC_COMM_SELF,k,k,NULL,&X0);
121:   MatDenseGetArray(X0,&pX0);
122:   for (i=0;i<k;i++) {
123:     PetscMemcpy(pX0+i*k,X+i*ldds,k*sizeof(PetscScalar));
124:   }
125:   MatDenseRestoreArray(X0,&pX0);
126:   BVSetActiveColumns(pep->V,0,k);
127:   BVMultInPlace(pep->V,X0,0,k);
128:   MatDestroy(&X0);
129:   BVSetActiveColumns(pep->V,0,k);
130:   DSRestoreArray(pep->ds,DS_MAT_X,&X);
131:   return(0);
132: }

136: /*
137:   Compute a step of Classical Gram-Schmidt orthogonalization
138: */
139: static PetscErrorCode PEPQArnoldiCGS(PEP pep,PetscScalar *H,PetscBLASInt ldh,PetscScalar *h,PetscBLASInt j,BV V,Vec t,Vec v,Vec w,PetscReal *onorm,PetscReal *norm,PetscScalar *work)
140: {
142:   PetscBLASInt   ione = 1,j_1 = j+1;
143:   PetscReal      x,y;
144:   PetscScalar    dot,one = 1.0,zero = 0.0;

147:   /* compute norm of v and w */
148:   if (onorm) {
149:     VecNorm(v,NORM_2,&x);
150:     VecNorm(w,NORM_2,&y);
151:     *onorm = PetscSqrtReal(x*x+y*y);
152:   }

154:   /* orthogonalize: compute h */
155:   BVDotVec(V,v,h);
156:   BVDotVec(V,w,work);
157:   if (j>0)
158:     PetscStackCallBLAS("BLASgemv",BLASgemv_("C",&j_1,&j,&one,H,&ldh,work,&ione,&one,h,&ione));
159:   VecDot(w,t,&dot);
160:   h[j] += dot;

162:   /* orthogonalize: update v and w */
163:   BVMultVec(V,-1.0,1.0,v,h);
164:   if (j>0) {
165:     PetscStackCallBLAS("BLASgemv",BLASgemv_("N",&j_1,&j,&one,H,&ldh,h,&ione,&zero,work,&ione));
166:     BVMultVec(V,-1.0,1.0,w,work);
167:   }
168:   VecAXPY(w,-h[j],t);

170:   /* compute norm of v and w */
171:   if (norm) {
172:     VecNorm(v,NORM_2,&x);
173:     VecNorm(w,NORM_2,&y);
174:     *norm = PetscSqrtReal(x*x+y*y);
175:   }
176:   return(0);
177: }

181: /*
182:   Compute a run of Q-Arnoldi iterations
183: */
184: static PetscErrorCode PEPQArnoldi(PEP pep,PetscScalar *H,PetscInt ldh,PetscInt k,PetscInt *M,Vec v,Vec w,PetscReal *beta,PetscBool *breakdown,PetscScalar *work)
185: {
186:   PetscErrorCode     ierr;
187:   PetscInt           i,j,l,m = *M;
188:   Vec                t = pep->work[2],u = pep->work[3];
189:   BVOrthogRefineType refinement;
190:   PetscReal          norm=0.0,onorm,eta;
191:   PetscScalar        *c = work + m;

194:   BVGetOrthogonalization(pep->V,NULL,&refinement,&eta,NULL);
195:   BVInsertVec(pep->V,k,v);
196:   for (j=k;j<m;j++) {
197:     /* apply operator */
198:     VecCopy(w,t);
199:     if (pep->Dr) {
200:       VecPointwiseMult(v,v,pep->Dr);
201:     }
202:     STMatMult(pep->st,0,v,u);
203:     VecCopy(t,v);
204:     if (pep->Dr) {
205:       VecPointwiseMult(t,t,pep->Dr);
206:     }
207:     STMatMult(pep->st,1,t,w);
208:     VecAXPY(u,pep->sfactor,w);
209:     STMatSolve(pep->st,u,w);
210:     VecScale(w,-1.0/(pep->sfactor*pep->sfactor));
211:     if (pep->Dr) {
212:       VecPointwiseDivide(w,w,pep->Dr);
213:     }
214:     VecCopy(v,t);
215:     BVSetActiveColumns(pep->V,0,j+1);

217:     /* orthogonalize */
218:     switch (refinement) {
219:       case BV_ORTHOG_REFINE_NEVER:
220:         PEPQArnoldiCGS(pep,H,ldh,H+ldh*j,j,pep->V,t,v,w,NULL,&norm,work);
221:         *breakdown = PETSC_FALSE;
222:         break;
223:       case BV_ORTHOG_REFINE_ALWAYS:
224:         PEPQArnoldiCGS(pep,H,ldh,H+ldh*j,j,pep->V,t,v,w,NULL,NULL,work);
225:         PEPQArnoldiCGS(pep,H,ldh,c,j,pep->V,t,v,w,&onorm,&norm,work);
226:         for (i=0;i<=j;i++) H[ldh*j+i] += c[i];
227:         if (norm < eta * onorm) *breakdown = PETSC_TRUE;
228:         else *breakdown = PETSC_FALSE;
229:         break;
230:       case BV_ORTHOG_REFINE_IFNEEDED:
231:         PEPQArnoldiCGS(pep,H,ldh,H+ldh*j,j,pep->V,t,v,w,&onorm,&norm,work);
232:         /* ||q|| < eta ||h|| */
233:         l = 1;
234:         while (l<3 && norm < eta * onorm) {
235:           l++;
236:           onorm = norm;
237:           PEPQArnoldiCGS(pep,H,ldh,c,j,pep->V,t,v,w,NULL,&norm,work);
238:           for (i=0;i<=j;i++) H[ldh*j+i] += c[i];
239:         }
240:         if (norm < eta * onorm) *breakdown = PETSC_TRUE;
241:         else *breakdown = PETSC_FALSE;
242:         break;
243:       default: SETERRQ(PetscObjectComm((PetscObject)pep),1,"Wrong value of ip->orth_ref");
244:     }
245:     VecScale(v,1.0/norm);
246:     VecScale(w,1.0/norm);

248:     H[j+1+ldh*j] = norm;
249:     if (j<m-1) {
250:       BVInsertVec(pep->V,j+1,v);
251:     }
252:   }
253:   *beta = norm;
254:   return(0);
255: }

259: PetscErrorCode PEPSolve_QArnoldi(PEP pep)
260: {
262:   PEP_QARNOLDI   *ctx = (PEP_QARNOLDI*)pep->data;
263:   PetscInt       j,k,l,lwork,nv,ld,newn,nconv;
264:   Vec            v=pep->work[0],w=pep->work[1];
265:   Mat            Q;
266:   PetscScalar    *S,*work;
267:   PetscReal      beta=0.0,norm,x,y;
268:   PetscBool      breakdown=PETSC_FALSE,sinv;

271:   DSGetLeadingDimension(pep->ds,&ld);
272:   lwork = 7*pep->ncv;
273:   PetscMalloc1(lwork,&work);
274:   PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv);
275:   RGPushScale(pep->rg,sinv?pep->sfactor:1.0/pep->sfactor);
276:   STScaleShift(pep->st,sinv?pep->sfactor:1.0/pep->sfactor);

278:   /* Get the starting Arnoldi vector */
279:   BVCopyVec(pep->V,0,v);
280:   BVCopyVec(pep->V,1,w);
281:   VecNorm(v,NORM_2,&x);
282:   VecNorm(w,NORM_2,&y);
283:   norm = PetscSqrtReal(x*x+y*y);
284:   VecScale(v,1.0/norm);
285:   VecScale(w,1.0/norm);

287:    /* Restart loop */
288:   l = 0;
289:   while (pep->reason == PEP_CONVERGED_ITERATING) {
290:     pep->its++;

292:     /* Compute an nv-step Arnoldi factorization */
293:     nv = PetscMin(pep->nconv+pep->mpd,pep->ncv);
294:     DSGetArray(pep->ds,DS_MAT_A,&S);
295:     PEPQArnoldi(pep,S,ld,pep->nconv+l,&nv,v,w,&beta,&breakdown,work);
296:     DSRestoreArray(pep->ds,DS_MAT_A,&S);
297:     DSSetDimensions(pep->ds,nv,0,pep->nconv,pep->nconv+l);
298:     if (l==0) {
299:       DSSetState(pep->ds,DS_STATE_INTERMEDIATE);
300:     } else {
301:       DSSetState(pep->ds,DS_STATE_RAW);
302:     }
303:     BVSetActiveColumns(pep->V,pep->nconv,nv);

305:     /* Solve projected problem */
306:     DSSolve(pep->ds,pep->eigr,pep->eigi);
307:     DSSort(pep->ds,pep->eigr,pep->eigi,NULL,NULL,NULL);
308:     DSUpdateExtraRow(pep->ds);

310:     /* Check convergence */
311:     PEPKrylovConvergence(pep,PETSC_FALSE,pep->nconv,nv-pep->nconv,beta,&k);
312:     (*pep->stopping)(pep,pep->its,pep->max_it,k,pep->nev,&pep->reason,pep->stoppingctx);
313:     nconv = k;

315:     /* Update l */
316:     if (pep->reason != PEP_CONVERGED_ITERATING || breakdown) l = 0;
317:     else l = PetscMax(1,(PetscInt)((nv-k)*ctx->keep));
318:     if (!ctx->lock && l>0) { l += k; k = 0; } /* non-locking variant: reset no. of converged pairs */

320:     if (pep->reason == PEP_CONVERGED_ITERATING) {
321:       if (breakdown) {
322:         /* Stop if breakdown */
323:         PetscInfo2(pep,"Breakdown Quadratic Arnoldi method (it=%D norm=%g)\n",pep->its,(double)beta);
324:         pep->reason = PEP_DIVERGED_BREAKDOWN;
325:       } else {
326:         /* Prepare the Rayleigh quotient for restart */
327:         DSTruncate(pep->ds,k+l);
328:         DSGetDimensions(pep->ds,&newn,NULL,NULL,NULL,NULL);
329:         l = newn-k;
330:       }
331:     }
332:     /* Update the corresponding vectors V(:,idx) = V*Q(:,idx) */
333:     DSGetMat(pep->ds,DS_MAT_Q,&Q);
334:     BVMultInPlace(pep->V,Q,pep->nconv,k+l);
335:     MatDestroy(&Q);

337:     pep->nconv = k;
338:     PEPMonitor(pep,pep->its,nconv,pep->eigr,pep->eigi,pep->errest,nv);
339:   }

341:   for (j=0;j<pep->nconv;j++) {
342:     pep->eigr[j] *= pep->sfactor;
343:     pep->eigi[j] *= pep->sfactor;
344:   }

346:   STScaleShift(pep->st,sinv?1.0/pep->sfactor:pep->sfactor);
347:   RGPopScale(pep->rg);

349:   /* truncate Schur decomposition and change the state to raw so that
350:      DSVectors() computes eigenvectors from scratch */
351:   DSSetDimensions(pep->ds,pep->nconv,0,0,0);
352:   DSSetState(pep->ds,DS_STATE_RAW);
353:   PetscFree(work);
354:   return(0);
355: }

359: static PetscErrorCode PEPQArnoldiSetRestart_QArnoldi(PEP pep,PetscReal keep)
360: {
361:   PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data;

364:   if (keep==PETSC_DEFAULT) ctx->keep = 0.5;
365:   else {
366:     if (keep<0.1 || keep>0.9) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"The keep argument must be in the range [0.1,0.9]");
367:     ctx->keep = keep;
368:   }
369:   return(0);
370: }

374: /*@
375:    PEPQArnoldiSetRestart - Sets the restart parameter for the Q-Arnoldi
376:    method, in particular the proportion of basis vectors that must be kept
377:    after restart.

379:    Logically Collective on PEP

381:    Input Parameters:
382: +  pep  - the eigenproblem solver context
383: -  keep - the number of vectors to be kept at restart

385:    Options Database Key:
386: .  -pep_qarnoldi_restart - Sets the restart parameter

388:    Notes:
389:    Allowed values are in the range [0.1,0.9]. The default is 0.5.

391:    Level: advanced

393: .seealso: PEPQArnoldiGetRestart()
394: @*/
395: PetscErrorCode PEPQArnoldiSetRestart(PEP pep,PetscReal keep)
396: {

402:   PetscTryMethod(pep,"PEPQArnoldiSetRestart_C",(PEP,PetscReal),(pep,keep));
403:   return(0);
404: }

408: static PetscErrorCode PEPQArnoldiGetRestart_QArnoldi(PEP pep,PetscReal *keep)
409: {
410:   PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data;

413:   *keep = ctx->keep;
414:   return(0);
415: }

419: /*@
420:    PEPQArnoldiGetRestart - Gets the restart parameter used in the Q-Arnoldi method.

422:    Not Collective

424:    Input Parameter:
425: .  pep - the eigenproblem solver context

427:    Output Parameter:
428: .  keep - the restart parameter

430:    Level: advanced

432: .seealso: PEPQArnoldiSetRestart()
433: @*/
434: PetscErrorCode PEPQArnoldiGetRestart(PEP pep,PetscReal *keep)
435: {

441:   PetscUseMethod(pep,"PEPQArnoldiGetRestart_C",(PEP,PetscReal*),(pep,keep));
442:   return(0);
443: }

447: static PetscErrorCode PEPQArnoldiSetLocking_QArnoldi(PEP pep,PetscBool lock)
448: {
449:   PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data;

452:   ctx->lock = lock;
453:   return(0);
454: }

458: /*@
459:    PEPQArnoldiSetLocking - Choose between locking and non-locking variants of
460:    the Q-Arnoldi method.

462:    Logically Collective on PEP

464:    Input Parameters:
465: +  pep  - the eigenproblem solver context
466: -  lock - true if the locking variant must be selected

468:    Options Database Key:
469: .  -pep_qarnoldi_locking - Sets the locking flag

471:    Notes:
472:    The default is to keep all directions in the working subspace even if
473:    already converged to working accuracy (the non-locking variant).
474:    This behaviour can be changed so that converged eigenpairs are locked
475:    when the method restarts.

477:    Note that the default behaviour is the opposite to Krylov solvers in EPS.

479:    Level: advanced

481: .seealso: PEPQArnoldiGetLocking()
482: @*/
483: PetscErrorCode PEPQArnoldiSetLocking(PEP pep,PetscBool lock)
484: {

490:   PetscTryMethod(pep,"PEPQArnoldiSetLocking_C",(PEP,PetscBool),(pep,lock));
491:   return(0);
492: }

496: static PetscErrorCode PEPQArnoldiGetLocking_QArnoldi(PEP pep,PetscBool *lock)
497: {
498:   PEP_QARNOLDI *ctx = (PEP_QARNOLDI*)pep->data;

501:   *lock = ctx->lock;
502:   return(0);
503: }

507: /*@
508:    PEPQArnoldiGetLocking - Gets the locking flag used in the Q-Arnoldi method.

510:    Not Collective

512:    Input Parameter:
513: .  pep - the eigenproblem solver context

515:    Output Parameter:
516: .  lock - the locking flag

518:    Level: advanced

520: .seealso: PEPQArnoldiSetLocking()
521: @*/
522: PetscErrorCode PEPQArnoldiGetLocking(PEP pep,PetscBool *lock)
523: {

529:   PetscUseMethod(pep,"PEPQArnoldiGetLocking_C",(PEP,PetscBool*),(pep,lock));
530:   return(0);
531: }

535: PetscErrorCode PEPSetFromOptions_QArnoldi(PetscOptionItems *PetscOptionsObject,PEP pep)
536: {
538:   PetscBool      flg,lock;
539:   PetscReal      keep;

542:   PetscOptionsHead(PetscOptionsObject,"PEP Q-Arnoldi Options");
543:   PetscOptionsReal("-pep_qarnoldi_restart","Proportion of vectors kept after restart","PEPQArnoldiSetRestart",0.5,&keep,&flg);
544:   if (flg) {
545:     PEPQArnoldiSetRestart(pep,keep);
546:   }
547:   PetscOptionsBool("-pep_qarnoldi_locking","Choose between locking and non-locking variants","PEPQArnoldiSetLocking",PETSC_FALSE,&lock,&flg);
548:   if (flg) {
549:     PEPQArnoldiSetLocking(pep,lock);
550:   }
551:   PetscOptionsTail();
552:   return(0);
553: }

557: PetscErrorCode PEPView_QArnoldi(PEP pep,PetscViewer viewer)
558: {
560:   PEP_QARNOLDI   *ctx = (PEP_QARNOLDI*)pep->data;
561:   PetscBool      isascii;

564:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
565:   if (isascii) {
566:     PetscViewerASCIIPrintf(viewer,"  Q-Arnoldi: %d%% of basis vectors kept after restart\n",(int)(100*ctx->keep));
567:     PetscViewerASCIIPrintf(viewer,"  Q-Arnoldi: using the %slocking variant\n",ctx->lock?"":"non-");
568:   }
569:   return(0);
570: }

574: PetscErrorCode PEPDestroy_QArnoldi(PEP pep)
575: {

579:   PetscFree(pep->data);
580:   PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiSetRestart_C",NULL);
581:   PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiGetRestart_C",NULL);
582:   PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiSetLocking_C",NULL);
583:   PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiGetLocking_C",NULL);
584:   return(0);
585: }

589: PETSC_EXTERN PetscErrorCode PEPCreate_QArnoldi(PEP pep)
590: {
591:   PEP_QARNOLDI   *ctx;

595:   PetscNewLog(pep,&ctx);
596:   pep->data = (void*)ctx;
597:   ctx->lock = PETSC_TRUE;

599:   pep->ops->solve          = PEPSolve_QArnoldi;
600:   pep->ops->setup          = PEPSetUp_QArnoldi;
601:   pep->ops->setfromoptions = PEPSetFromOptions_QArnoldi;
602:   pep->ops->destroy        = PEPDestroy_QArnoldi;
603:   pep->ops->view           = PEPView_QArnoldi;
604:   pep->ops->backtransform  = PEPBackTransform_Default;
605:   pep->ops->computevectors = PEPComputeVectors_Default;
606:   pep->ops->computevectors = PEPExtractVectors_QArnoldi;
607:   PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiSetRestart_C",PEPQArnoldiSetRestart_QArnoldi);
608:   PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiGetRestart_C",PEPQArnoldiGetRestart_QArnoldi);
609:   PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiSetLocking_C",PEPQArnoldiSetLocking_QArnoldi);
610:   PetscObjectComposeFunction((PetscObject)pep,"PEPQArnoldiGetLocking_C",PEPQArnoldiGetLocking_QArnoldi);
611:   return(0);
612: }