Actual source code: jd.c

slepc-3.18.0 2022-10-01
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, 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 eigensolver: "jd"

 13:    Method: Jacobi-Davidson

 15:    Algorithm:

 17:        Jacobi-Davidson with various subspace extraction and
 18:        restart techniques.

 20:    References:

 22:        [1] G.L.G. Sleijpen and H.A. van der Vorst, "A Jacobi-Davidson
 23:            iteration method for linear eigenvalue problems", SIAM J.
 24:            Matrix Anal. Appl. 17(2):401-425, 1996.

 26:        [2] E. Romero and J.E. Roman, "A parallel implementation of
 27:            Davidson methods for large-scale eigenvalue problems in
 28:            SLEPc", ACM Trans. Math. Software 40(2), Article 13, 2014.
 29: */

 31: #include <slepc/private/epsimpl.h>
 32: #include <../src/eps/impls/davidson/davidson.h>

 34: PetscErrorCode EPSSetFromOptions_JD(EPS eps,PetscOptionItems *PetscOptionsObject)
 35: {
 36:   PetscBool      flg,flg2,op,orth;
 37:   PetscInt       opi,opi0;
 38:   PetscReal      opf;

 40:   PetscOptionsHeadBegin(PetscOptionsObject,"EPS Jacobi-Davidson (JD) Options");

 42:     EPSJDGetKrylovStart(eps,&op);
 43:     PetscOptionsBool("-eps_jd_krylov_start","Start the search subspace with a Krylov basis","EPSJDSetKrylovStart",op,&op,&flg);
 44:     if (flg) EPSJDSetKrylovStart(eps,op);

 46:     EPSJDGetBOrth(eps,&orth);
 47:     PetscOptionsBool("-eps_jd_borth","Use B-orthogonalization in the search subspace","EPSJDSetBOrth",op,&op,&flg);
 48:     if (flg) EPSJDSetBOrth(eps,op);

 50:     EPSJDGetBlockSize(eps,&opi);
 51:     PetscOptionsInt("-eps_jd_blocksize","Number of vectors to add to the search subspace","EPSJDSetBlockSize",opi,&opi,&flg);
 52:     if (flg) EPSJDSetBlockSize(eps,opi);

 54:     EPSJDGetRestart(eps,&opi,&opi0);
 55:     PetscOptionsInt("-eps_jd_minv","Size of the search subspace after restarting","EPSJDSetRestart",opi,&opi,&flg);
 56:     PetscOptionsInt("-eps_jd_plusk","Number of eigenvectors saved from the previous iteration when restarting","EPSJDSetRestart",opi0,&opi0,&flg2);
 57:     if (flg || flg2) EPSJDSetRestart(eps,opi,opi0);

 59:     EPSJDGetInitialSize(eps,&opi);
 60:     PetscOptionsInt("-eps_jd_initial_size","Initial size of the search subspace","EPSJDSetInitialSize",opi,&opi,&flg);
 61:     if (flg) EPSJDSetInitialSize(eps,opi);

 63:     EPSJDGetFix(eps,&opf);
 64:     PetscOptionsReal("-eps_jd_fix","Tolerance for changing the target in the correction equation","EPSJDSetFix",opf,&opf,&flg);
 65:     if (flg) EPSJDSetFix(eps,opf);

 67:     EPSJDGetConstCorrectionTol(eps,&op);
 68:     PetscOptionsBool("-eps_jd_const_correction_tol","Disable the dynamic stopping criterion when solving the correction equation","EPSJDSetConstCorrectionTol",op,&op,&flg);
 69:     if (flg) EPSJDSetConstCorrectionTol(eps,op);

 71:   PetscOptionsHeadEnd();
 72:   return 0;
 73: }

 75: PetscErrorCode EPSSetDefaultST_JD(EPS eps)
 76: {
 77:   KSP            ksp;

 79:   if (!((PetscObject)eps->st)->type_name) {
 80:     STSetType(eps->st,STPRECOND);
 81:     STPrecondSetKSPHasMat(eps->st,PETSC_TRUE);
 82:   }
 83:   STGetKSP(eps->st,&ksp);
 84:   if (!((PetscObject)ksp)->type_name) {
 85:     KSPSetType(ksp,KSPBCGSL);
 86:     KSPSetTolerances(ksp,1e-4,PETSC_DEFAULT,PETSC_DEFAULT,90);
 87:   }
 88:   return 0;
 89: }

 91: PetscErrorCode EPSSetUp_JD(EPS eps)
 92: {
 93:   PetscBool      t;
 94:   KSP            ksp;

 96:   /* Setup common for all davidson solvers */
 97:   EPSSetUp_XD(eps);

 99:   /* Check some constraints */
100:   STGetKSP(eps->st,&ksp);
101:   PetscObjectTypeCompare((PetscObject)ksp,KSPPREONLY,&t);
103:   return 0;
104: }

106: PetscErrorCode EPSView_JD(EPS eps,PetscViewer viewer)
107: {
108:   PetscBool      isascii,opb;
109:   PetscReal      opf;
110:   PetscInt       opi,opi0;

112:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
113:   if (isascii) {
114:     EPSXDGetBOrth_XD(eps,&opb);
115:     if (opb) PetscViewerASCIIPrintf(viewer,"  search subspace is B-orthogonalized\n");
116:     else PetscViewerASCIIPrintf(viewer,"  search subspace is orthogonalized\n");
117:     EPSXDGetBlockSize_XD(eps,&opi);
118:     PetscViewerASCIIPrintf(viewer,"  block size=%" PetscInt_FMT "\n",opi);
119:     EPSXDGetKrylovStart_XD(eps,&opb);
120:     if (!opb) PetscViewerASCIIPrintf(viewer,"  type of the initial subspace: non-Krylov\n");
121:     else PetscViewerASCIIPrintf(viewer,"  type of the initial subspace: Krylov\n");
122:     EPSXDGetRestart_XD(eps,&opi,&opi0);
123:     PetscViewerASCIIPrintf(viewer,"  size of the subspace after restarting: %" PetscInt_FMT "\n",opi);
124:     PetscViewerASCIIPrintf(viewer,"  number of vectors after restarting from the previous iteration: %" PetscInt_FMT "\n",opi0);

126:     EPSJDGetFix_JD(eps,&opf);
127:     PetscViewerASCIIPrintf(viewer,"  threshold for changing the target in the correction equation (fix): %g\n",(double)opf);

129:     EPSJDGetConstCorrectionTol_JD(eps,&opb);
130:     if (!opb) PetscViewerASCIIPrintf(viewer,"  using dynamic tolerance for the correction equation\n");
131:   }
132:   return 0;
133: }

135: PetscErrorCode EPSDestroy_JD(EPS eps)
136: {
137:   PetscFree(eps->data);
138:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",NULL);
139:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",NULL);
140:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",NULL);
141:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",NULL);
142:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",NULL);
143:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",NULL);
144:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",NULL);
145:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",NULL);
146:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",NULL);
147:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",NULL);
148:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",NULL);
149:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",NULL);
150:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",NULL);
151:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",NULL);
152:   return 0;
153: }

155: /*@
156:    EPSJDSetKrylovStart - Activates or deactivates starting the searching
157:    subspace with a Krylov basis.

159:    Logically Collective on eps

161:    Input Parameters:
162: +  eps - the eigenproblem solver context
163: -  krylovstart - boolean flag

165:    Options Database Key:
166: .  -eps_jd_krylov_start - Activates starting the searching subspace with a
167:     Krylov basis

169:    Level: advanced

171: .seealso: EPSJDGetKrylovStart()
172: @*/
173: PetscErrorCode EPSJDSetKrylovStart(EPS eps,PetscBool krylovstart)
174: {
177:   PetscTryMethod(eps,"EPSJDSetKrylovStart_C",(EPS,PetscBool),(eps,krylovstart));
178:   return 0;
179: }

181: /*@
182:    EPSJDGetKrylovStart - Returns a flag indicating if the searching subspace is started with a
183:    Krylov basis.

185:    Not Collective

187:    Input Parameter:
188: .  eps - the eigenproblem solver context

190:    Output Parameters:
191: .  krylovstart - boolean flag indicating if the searching subspace is started
192:    with a Krylov basis

194:    Level: advanced

196: .seealso: EPSJDSetKrylovStart()
197: @*/
198: PetscErrorCode EPSJDGetKrylovStart(EPS eps,PetscBool *krylovstart)
199: {
202:   PetscUseMethod(eps,"EPSJDGetKrylovStart_C",(EPS,PetscBool*),(eps,krylovstart));
203:   return 0;
204: }

206: /*@
207:    EPSJDSetBlockSize - Sets the number of vectors to be added to the searching space
208:    in every iteration.

210:    Logically Collective on eps

212:    Input Parameters:
213: +  eps - the eigenproblem solver context
214: -  blocksize - number of vectors added to the search space in every iteration

216:    Options Database Key:
217: .  -eps_jd_blocksize - number of vectors added to the searching space every iteration

219:    Level: advanced

221: .seealso: EPSJDSetKrylovStart()
222: @*/
223: PetscErrorCode EPSJDSetBlockSize(EPS eps,PetscInt blocksize)
224: {
227:   PetscTryMethod(eps,"EPSJDSetBlockSize_C",(EPS,PetscInt),(eps,blocksize));
228:   return 0;
229: }

231: /*@
232:    EPSJDGetBlockSize - Returns the number of vectors to be added to the searching space
233:    in every iteration.

235:    Not Collective

237:    Input Parameter:
238: .  eps - the eigenproblem solver context

240:    Output Parameter:
241: .  blocksize - number of vectors added to the search space in every iteration

243:    Level: advanced

245: .seealso: EPSJDSetBlockSize()
246: @*/
247: PetscErrorCode EPSJDGetBlockSize(EPS eps,PetscInt *blocksize)
248: {
251:   PetscUseMethod(eps,"EPSJDGetBlockSize_C",(EPS,PetscInt*),(eps,blocksize));
252:   return 0;
253: }

255: /*@
256:    EPSJDSetRestart - Sets the number of vectors of the searching space after
257:    restarting and the number of vectors saved from the previous iteration.

259:    Logically Collective on eps

261:    Input Parameters:
262: +  eps - the eigenproblem solver context
263: .  minv - number of vectors of the searching subspace after restarting
264: -  plusk - number of vectors saved from the previous iteration

266:    Options Database Keys:
267: +  -eps_jd_minv - number of vectors of the searching subspace after restarting
268: -  -eps_jd_plusk - number of vectors saved from the previous iteration

270:    Level: advanced

272: .seealso: EPSJDGetRestart()
273: @*/
274: PetscErrorCode EPSJDSetRestart(EPS eps,PetscInt minv,PetscInt plusk)
275: {
279:   PetscTryMethod(eps,"EPSJDSetRestart_C",(EPS,PetscInt,PetscInt),(eps,minv,plusk));
280:   return 0;
281: }

283: /*@
284:    EPSJDGetRestart - Gets the number of vectors of the searching space after
285:    restarting and the number of vectors saved from the previous iteration.

287:    Not Collective

289:    Input Parameter:
290: .  eps - the eigenproblem solver context

292:    Output Parameters:
293: +  minv - number of vectors of the searching subspace after restarting
294: -  plusk - number of vectors saved from the previous iteration

296:    Level: advanced

298: .seealso: EPSJDSetRestart()
299: @*/
300: PetscErrorCode EPSJDGetRestart(EPS eps,PetscInt *minv,PetscInt *plusk)
301: {
303:   PetscUseMethod(eps,"EPSJDGetRestart_C",(EPS,PetscInt*,PetscInt*),(eps,minv,plusk));
304:   return 0;
305: }

307: /*@
308:    EPSJDSetInitialSize - Sets the initial size of the searching space.

310:    Logically Collective on eps

312:    Input Parameters:
313: +  eps - the eigenproblem solver context
314: -  initialsize - number of vectors of the initial searching subspace

316:    Options Database Key:
317: .  -eps_jd_initial_size - number of vectors of the initial searching subspace

319:    Notes:
320:    If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
321:    EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
322:    provided vectors are not enough, the solver completes the subspace with
323:    random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
324:    gets the first vector provided by the user or, if not available, a random vector,
325:    and expands the Krylov basis up to initialsize vectors.

327:    Level: advanced

329: .seealso: EPSJDGetInitialSize(), EPSJDGetKrylovStart()
330: @*/
331: PetscErrorCode EPSJDSetInitialSize(EPS eps,PetscInt initialsize)
332: {
335:   PetscTryMethod(eps,"EPSJDSetInitialSize_C",(EPS,PetscInt),(eps,initialsize));
336:   return 0;
337: }

339: /*@
340:    EPSJDGetInitialSize - Returns the initial size of the searching space.

342:    Not Collective

344:    Input Parameter:
345: .  eps - the eigenproblem solver context

347:    Output Parameter:
348: .  initialsize - number of vectors of the initial searching subspace

350:    Notes:
351:    If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
352:    EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
353:    provided vectors are not enough, the solver completes the subspace with
354:    random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
355:    gets the first vector provided by the user or, if not available, a random vector,
356:    and expands the Krylov basis up to initialsize vectors.

358:    Level: advanced

360: .seealso: EPSJDSetInitialSize(), EPSJDGetKrylovStart()
361: @*/
362: PetscErrorCode EPSJDGetInitialSize(EPS eps,PetscInt *initialsize)
363: {
366:   PetscUseMethod(eps,"EPSJDGetInitialSize_C",(EPS,PetscInt*),(eps,initialsize));
367:   return 0;
368: }

370: PetscErrorCode EPSJDSetFix_JD(EPS eps,PetscReal fix)
371: {
372:   EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;

374:   if (fix == PETSC_DEFAULT || fix == PETSC_DECIDE) fix = 0.01;
376:   data->fix = fix;
377:   return 0;
378: }

380: /*@
381:    EPSJDSetFix - Sets the threshold for changing the target in the correction
382:    equation.

384:    Logically Collective on eps

386:    Input Parameters:
387: +  eps - the eigenproblem solver context
388: -  fix - threshold for changing the target

390:    Options Database Key:
391: .  -eps_jd_fix - the fix value

393:    Note:
394:    The target in the correction equation is fixed at the first iterations.
395:    When the norm of the residual vector is lower than the fix value,
396:    the target is set to the corresponding eigenvalue.

398:    Level: advanced

400: .seealso: EPSJDGetFix()
401: @*/
402: PetscErrorCode EPSJDSetFix(EPS eps,PetscReal fix)
403: {
406:   PetscTryMethod(eps,"EPSJDSetFix_C",(EPS,PetscReal),(eps,fix));
407:   return 0;
408: }

410: PetscErrorCode EPSJDGetFix_JD(EPS eps,PetscReal *fix)
411: {
412:   EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;

414:   *fix = data->fix;
415:   return 0;
416: }

418: /*@
419:    EPSJDGetFix - Returns the threshold for changing the target in the correction
420:    equation.

422:    Not Collective

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

427:    Output Parameter:
428: .  fix - threshold for changing the target

430:    Note:
431:    The target in the correction equation is fixed at the first iterations.
432:    When the norm of the residual vector is lower than the fix value,
433:    the target is set to the corresponding eigenvalue.

435:    Level: advanced

437: .seealso: EPSJDSetFix()
438: @*/
439: PetscErrorCode EPSJDGetFix(EPS eps,PetscReal *fix)
440: {
443:   PetscUseMethod(eps,"EPSJDGetFix_C",(EPS,PetscReal*),(eps,fix));
444:   return 0;
445: }

447: PetscErrorCode EPSJDSetConstCorrectionTol_JD(EPS eps,PetscBool constant)
448: {
449:   EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;

451:   data->dynamic = PetscNot(constant);
452:   return 0;
453: }

455: /*@
456:    EPSJDSetConstCorrectionTol - If true, deactivates the dynamic stopping criterion
457:    (also called Newton) that sets the KSP relative tolerance
458:    to 0.5**i, where i is the number of EPS iterations from the last converged value.

460:    Logically Collective on eps

462:    Input Parameters:
463: +  eps - the eigenproblem solver context
464: -  constant - if false, the KSP relative tolerance is set to 0.5**i.

466:    Options Database Key:
467: .  -eps_jd_const_correction_tol - Deactivates the dynamic stopping criterion.

469:    Level: advanced

471: .seealso: EPSJDGetConstCorrectionTol()
472: @*/
473: PetscErrorCode EPSJDSetConstCorrectionTol(EPS eps,PetscBool constant)
474: {
477:   PetscTryMethod(eps,"EPSJDSetConstCorrectionTol_C",(EPS,PetscBool),(eps,constant));
478:   return 0;
479: }

481: PetscErrorCode EPSJDGetConstCorrectionTol_JD(EPS eps,PetscBool *constant)
482: {
483:   EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;

485:   *constant = PetscNot(data->dynamic);
486:   return 0;
487: }

489: /*@
490:    EPSJDGetConstCorrectionTol - Returns a flag indicating if the dynamic stopping is being used for
491:    solving the correction equation.

493:    Not Collective

495:    Input Parameter:
496: .  eps - the eigenproblem solver context

498:    Output Parameter:
499: .  constant - boolean flag indicating if the dynamic stopping criterion is not being used.

501:    Notes:
502:    If the flag is false the KSP relative tolerance is set to 0.5**i, where i is the number
503:    of EPS iterations from the last converged value.

505:    Level: advanced

507: .seealso: EPSJDSetConstCorrectionTol()
508: @*/
509: PetscErrorCode EPSJDGetConstCorrectionTol(EPS eps,PetscBool *constant)
510: {
513:   PetscUseMethod(eps,"EPSJDGetConstCorrectionTol_C",(EPS,PetscBool*),(eps,constant));
514:   return 0;
515: }

517: /*@
518:    EPSJDSetBOrth - Selects the orthogonalization that will be used in the search
519:    subspace in case of generalized Hermitian problems.

521:    Logically Collective on eps

523:    Input Parameters:
524: +  eps   - the eigenproblem solver context
525: -  borth - whether to B-orthogonalize the search subspace

527:    Options Database Key:
528: .  -eps_jd_borth - Set the orthogonalization used in the search subspace

530:    Level: advanced

532: .seealso: EPSJDGetBOrth()
533: @*/
534: PetscErrorCode EPSJDSetBOrth(EPS eps,PetscBool borth)
535: {
538:   PetscTryMethod(eps,"EPSJDSetBOrth_C",(EPS,PetscBool),(eps,borth));
539:   return 0;
540: }

542: /*@
543:    EPSJDGetBOrth - Returns the orthogonalization used in the search
544:    subspace in case of generalized Hermitian problems.

546:    Not Collective

548:    Input Parameter:
549: .  eps - the eigenproblem solver context

551:    Output Parameters:
552: .  borth - whether to B-orthogonalize the search subspace

554:    Level: advanced

556: .seealso: EPSJDSetBOrth()
557: @*/
558: PetscErrorCode EPSJDGetBOrth(EPS eps,PetscBool *borth)
559: {
562:   PetscUseMethod(eps,"EPSJDGetBOrth_C",(EPS,PetscBool*),(eps,borth));
563:   return 0;
564: }

566: SLEPC_EXTERN PetscErrorCode EPSCreate_JD(EPS eps)
567: {
568:   EPS_DAVIDSON   *data;

570:   PetscNew(&data);
571:   eps->data = (void*)data;

573:   data->blocksize   = 1;
574:   data->initialsize = 0;
575:   data->minv        = 0;
576:   data->plusk       = PETSC_DEFAULT;
577:   data->ipB         = PETSC_TRUE;
578:   data->fix         = 0.01;
579:   data->krylovstart = PETSC_FALSE;
580:   data->dynamic     = PETSC_FALSE;

582:   eps->useds = PETSC_TRUE;
583:   eps->categ = EPS_CATEGORY_PRECOND;

585:   eps->ops->solve          = EPSSolve_XD;
586:   eps->ops->setup          = EPSSetUp_JD;
587:   eps->ops->setupsort      = EPSSetUpSort_Default;
588:   eps->ops->setfromoptions = EPSSetFromOptions_JD;
589:   eps->ops->destroy        = EPSDestroy_JD;
590:   eps->ops->reset          = EPSReset_XD;
591:   eps->ops->view           = EPSView_JD;
592:   eps->ops->backtransform  = EPSBackTransform_Default;
593:   eps->ops->computevectors = EPSComputeVectors_XD;
594:   eps->ops->setdefaultst   = EPSSetDefaultST_JD;

596:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",EPSXDSetKrylovStart_XD);
597:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",EPSXDGetKrylovStart_XD);
598:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",EPSXDSetBlockSize_XD);
599:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",EPSXDGetBlockSize_XD);
600:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",EPSXDSetRestart_XD);
601:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",EPSXDGetRestart_XD);
602:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",EPSXDSetInitialSize_XD);
603:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",EPSXDGetInitialSize_XD);
604:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",EPSJDSetFix_JD);
605:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",EPSJDGetFix_JD);
606:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",EPSJDSetConstCorrectionTol_JD);
607:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",EPSJDGetConstCorrectionTol_JD);
608:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",EPSXDSetBOrth_XD);
609:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",EPSXDGetBOrth_XD);
610:   return 0;
611: }