Actual source code: jd.c

slepc-3.7.0 2016-05-16
Report Typos and Errors
  1: /*

  3:    SLEPc eigensolver: "jd"

  5:    Method: Jacobi-Davidson

  7:    Algorithm:

  9:        Jacobi-Davidson with various subspace extraction and
 10:        restart techniques.

 12:    References:

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

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

 22:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 23:    SLEPc - Scalable Library for Eigenvalue Problem Computations
 24:    Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain

 26:    This file is part of SLEPc.

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

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

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

 42: #include <slepc/private/epsimpl.h>                /*I "slepceps.h" I*/
 43: #include <../src/eps/impls/davidson/davidson.h>

 47: PetscErrorCode EPSSetFromOptions_JD(PetscOptionItems *PetscOptionsObject,EPS eps)
 48: {
 50:   PetscBool      flg,op;
 51:   PetscInt       opi,opi0;
 52:   PetscReal      opf;
 53:   KSP            ksp;
 54:   PetscBool      orth;
 55:   const char     *orth_list[2] = {"I","B"};

 58:   PetscOptionsHead(PetscOptionsObject,"EPS Jacobi-Davidson (JD) Options");

 60:   EPSJDGetKrylovStart(eps,&op);
 61:   PetscOptionsBool("-eps_jd_krylov_start","Start the searching subspace with a krylov basis","EPSJDSetKrylovStart",op,&op,&flg);
 62:   if (flg) { EPSJDSetKrylovStart(eps,op); }

 64:   EPSJDGetBlockSize(eps,&opi);
 65:   PetscOptionsInt("-eps_jd_blocksize","Number vectors add to the searching subspace","EPSJDSetBlockSize",opi,&opi,&flg);
 66:   if (flg) { EPSJDSetBlockSize(eps,opi); }

 68:   EPSJDGetRestart(eps,&opi,&opi0);
 69:   PetscOptionsInt("-eps_jd_minv","Set the size of the searching subspace after restarting","EPSJDSetRestart",opi,&opi,&flg);
 70:   if (flg) { EPSJDSetRestart(eps,opi,opi0); }

 72:   PetscOptionsInt("-eps_jd_plusk","Set the number of saved eigenvectors from the previous iteration when restarting","EPSJDSetRestart",opi0,&opi0,&flg);
 73:   if (flg) { EPSJDSetRestart(eps,opi,opi0); }

 75:   EPSJDGetInitialSize(eps,&opi);
 76:   PetscOptionsInt("-eps_jd_initial_size","Set the initial size of the searching subspace","EPSJDSetInitialSize",opi,&opi,&flg);
 77:   if (flg) { EPSJDSetInitialSize(eps,opi); }

 79:   EPSJDGetFix(eps,&opf);
 80:   PetscOptionsReal("-eps_jd_fix","Set the tolerance for changing the target in the correction equation","EPSJDSetFix",opf,&opf,&flg);
 81:   if (flg) { EPSJDSetFix(eps,opf); }

 83:   EPSJDGetBOrth(eps,&orth);
 84:   PetscOptionsEList("-eps_jd_borth","orthogonalization used in the search subspace","EPSJDSetBOrth",orth_list,2,orth_list[orth?1:0],&opi,&flg);
 85:   if (flg) { EPSJDSetBOrth(eps,opi==1?PETSC_TRUE:PETSC_FALSE); }

 87:   EPSJDGetConstCorrectionTol(eps,&op);
 88:   PetscOptionsBool("-eps_jd_const_correction_tol","Disable the dynamic stopping criterion when solving the correction equation","EPSJDSetConstCorrectionTol",op,&op,&flg);
 89:   if (flg) { EPSJDSetConstCorrectionTol(eps,op); }

 91:   EPSJDGetWindowSizes(eps,&opi,&opi0);
 92:   PetscOptionsInt("-eps_jd_pwindow","(Experimental!) Set the number of converged vectors in the projector","EPSJDSetWindowSizes",opi,&opi,&flg);
 93:   if (flg) { EPSJDSetWindowSizes(eps,opi,opi0); }

 95:   PetscOptionsInt("-eps_jd_qwindow","(Experimental!) Set the number of converged vectors in the projected problem","EPSJDSetWindowSizes",opi0,&opi0,&flg);
 96:   if (flg) { EPSJDSetWindowSizes(eps,opi,opi0); }

 98:   /* Set STPrecond as the default ST */
 99:   if (!((PetscObject)eps->st)->type_name) {
100:     STSetType(eps->st,STPRECOND);
101:   }
102:   STPrecondSetKSPHasMat(eps->st,PETSC_FALSE);

104:   /* Set the default options of the KSP */
105:   STGetKSP(eps->st,&ksp);
106:   if (!((PetscObject)ksp)->type_name) {
107:     KSPSetType(ksp,KSPBCGSL);
108:     KSPSetTolerances(ksp,1e-4,PETSC_DEFAULT,PETSC_DEFAULT,90);
109:   }
110:   PetscOptionsTail();
111:   return(0);
112: }

116: PetscErrorCode EPSSetUp_JD(EPS eps)
117: {
119:   PetscBool      t;
120:   KSP            ksp;

123:   /* Setup common for all davidson solvers */
124:   EPSSetUp_XD(eps);

126:   /* Set the default options of the KSP */
127:   STGetKSP(eps->st,&ksp);
128:   if (!((PetscObject)ksp)->type_name) {
129:     KSPSetType(ksp,KSPBCGSL);
130:     KSPSetTolerances(ksp,1e-4,PETSC_DEFAULT,PETSC_DEFAULT,90);
131:   }

133:   /* Check some constraints */
134:   PetscObjectTypeCompare((PetscObject)ksp,KSPPREONLY,&t);
135:   if (t) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"EPSJD does not work with KSPPREONLY");
136:   return(0);
137: }

141: PetscErrorCode EPSView_JD(EPS eps,PetscViewer viewer)
142: {
144:   PetscBool      isascii,opb;
145:   PetscInt       opi,opi0;
146:   PetscBool      borth;

149:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
150:   if (isascii) {
151:     EPSXDGetBOrth_XD(eps,&borth);
152:     if (borth) {
153:       PetscViewerASCIIPrintf(viewer,"  JD: search subspace is B-orthogonalized\n");
154:     } else {
155:       PetscViewerASCIIPrintf(viewer,"  JD: search subspace is orthogonalized\n");
156:     }
157:     EPSXDGetBlockSize_XD(eps,&opi);
158:     PetscViewerASCIIPrintf(viewer,"  JD: block size=%D\n",opi);
159:     EPSXDGetKrylovStart_XD(eps,&opb);
160:     if (!opb) {
161:       PetscViewerASCIIPrintf(viewer,"  JD: type of the initial subspace: non-Krylov\n");
162:     } else {
163:       PetscViewerASCIIPrintf(viewer,"  JD: type of the initial subspace: Krylov\n");
164:     }
165:     EPSXDGetRestart_XD(eps,&opi,&opi0);
166:     PetscViewerASCIIPrintf(viewer,"  JD: size of the subspace after restarting: %D\n",opi);
167:     PetscViewerASCIIPrintf(viewer,"  JD: number of vectors after restarting from the previous iteration: %D\n",opi0);
168:   }
169:   return(0);
170: }

174: PetscErrorCode EPSDestroy_JD(EPS eps)
175: {
176:   PetscErrorCode  ierr;

179:   PetscFree(eps->data);
180:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",NULL);
181:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",NULL);
182:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",NULL);
183:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",NULL);
184:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",NULL);
185:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",NULL);
186:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",NULL);
187:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",NULL);
188:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",NULL);
189:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",NULL);
190:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",NULL);
191:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",NULL);
192:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetWindowSizes_C",NULL);
193:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetWindowSizes_C",NULL);
194:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",NULL);
195:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",NULL);
196:   return(0);
197: }

201: /*@
202:    EPSJDSetKrylovStart - Activates or deactivates starting the searching
203:    subspace with a Krylov basis.

205:    Logically Collective on EPS

207:    Input Parameters:
208: +  eps - the eigenproblem solver context
209: -  krylovstart - boolean flag

211:    Options Database Key:
212: .  -eps_jd_krylov_start - Activates starting the searching subspace with a
213:     Krylov basis

215:    Level: advanced

217: .seealso: EPSJDGetKrylovStart()
218: @*/
219: PetscErrorCode EPSJDSetKrylovStart(EPS eps,PetscBool krylovstart)
220: {

226:   PetscTryMethod(eps,"EPSJDSetKrylovStart_C",(EPS,PetscBool),(eps,krylovstart));
227:   return(0);
228: }

232: /*@
233:    EPSJDGetKrylovStart - Returns a flag indicating if the searching subspace is started with a
234:    Krylov basis.

236:    Not Collective

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

241:    Output Parameters:
242: .  krylovstart - boolean flag indicating if the searching subspace is started
243:    with a Krylov basis

245:    Level: advanced

247: .seealso: EPSJDGetKrylovStart()
248: @*/
249: PetscErrorCode EPSJDGetKrylovStart(EPS eps,PetscBool *krylovstart)
250: {

256:   PetscUseMethod(eps,"EPSJDGetKrylovStart_C",(EPS,PetscBool*),(eps,krylovstart));
257:   return(0);
258: }

262: /*@
263:    EPSJDSetBlockSize - Sets the number of vectors to be added to the searching space
264:    in every iteration.

266:    Logically Collective on EPS

268:    Input Parameters:
269: +  eps - the eigenproblem solver context
270: -  blocksize - number of vectors added to the search space in every iteration

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

275:    Level: advanced

277: .seealso: EPSJDSetKrylovStart()
278: @*/
279: PetscErrorCode EPSJDSetBlockSize(EPS eps,PetscInt blocksize)
280: {

286:   PetscTryMethod(eps,"EPSJDSetBlockSize_C",(EPS,PetscInt),(eps,blocksize));
287:   return(0);
288: }

292: /*@
293:    EPSJDGetBlockSize - Returns the number of vectors to be added to the searching space
294:    in every iteration.

296:    Not Collective

298:    Input Parameter:
299: .  eps - the eigenproblem solver context

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

304:    Level: advanced

306: .seealso: EPSJDSetBlockSize()
307: @*/
308: PetscErrorCode EPSJDGetBlockSize(EPS eps,PetscInt *blocksize)
309: {

315:   PetscUseMethod(eps,"EPSJDGetBlockSize_C",(EPS,PetscInt*),(eps,blocksize));
316:   return(0);
317: }

321: /*@
322:    EPSJDGetRestart - Gets the number of vectors of the searching space after
323:    restarting and the number of vectors saved from the previous iteration.

325:    Not Collective

327:    Input Parameter:
328: .  eps - the eigenproblem solver context

330:    Output Parameter:
331: +  minv - number of vectors of the searching subspace after restarting
332: -  plusk - number of vectors saved from the previous iteration

334:    Level: advanced

336: .seealso: EPSJDSetRestart()
337: @*/
338: PetscErrorCode EPSJDGetRestart(EPS eps,PetscInt *minv,PetscInt *plusk)
339: {

344:   PetscUseMethod(eps,"EPSJDGetRestart_C",(EPS,PetscInt*,PetscInt*),(eps,minv,plusk));
345:   return(0);
346: }

350: /*@
351:    EPSJDSetRestart - Sets the number of vectors of the searching space after
352:    restarting and the number of vectors saved from the previous iteration.

354:    Logically Collective on EPS

356:    Input Parameters:
357: +  eps - the eigenproblem solver context
358: .  minv - number of vectors of the searching subspace after restarting
359: -  plusk - number of vectors saved from the previous iteration

361:    Options Database Keys:
362: +  -eps_jd_minv - number of vectors of the searching subspace after restarting
363: -  -eps_jd_plusk - number of vectors saved from the previous iteration

365:    Level: advanced

367: .seealso: EPSJDGetRestart()
368: @*/
369: PetscErrorCode EPSJDSetRestart(EPS eps,PetscInt minv,PetscInt plusk)
370: {

377:   PetscTryMethod(eps,"EPSJDSetRestart_C",(EPS,PetscInt,PetscInt),(eps,minv,plusk));
378:   return(0);
379: }

383: /*@
384:    EPSJDGetInitialSize - Returns the initial size of the searching space.

386:    Not Collective

388:    Input Parameter:
389: .  eps - the eigenproblem solver context

391:    Output Parameter:
392: .  initialsize - number of vectors of the initial searching subspace

394:    Notes:
395:    If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
396:    EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
397:    provided vectors are not enough, the solver completes the subspace with
398:    random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
399:    gets the first vector provided by the user or, if not available, a random vector,
400:    and expands the Krylov basis up to initialsize vectors.

402:    Level: advanced

404: .seealso: EPSJDSetInitialSize(), EPSJDGetKrylovStart()
405: @*/
406: PetscErrorCode EPSJDGetInitialSize(EPS eps,PetscInt *initialsize)
407: {

413:   PetscUseMethod(eps,"EPSJDGetInitialSize_C",(EPS,PetscInt*),(eps,initialsize));
414:   return(0);
415: }

419: /*@
420:    EPSJDSetInitialSize - Sets the initial size of the searching space.

422:    Logically Collective on EPS

424:    Input Parameters:
425: +  eps - the eigenproblem solver context
426: -  initialsize - number of vectors of the initial searching subspace

428:    Options Database Key:
429: .  -eps_jd_initial_size - number of vectors of the initial searching subspace

431:    Notes:
432:    If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
433:    EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
434:    provided vectors are not enough, the solver completes the subspace with
435:    random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
436:    gets the first vector provided by the user or, if not available, a random vector,
437:    and expands the Krylov basis up to initialsize vectors.

439:    Level: advanced

441: .seealso: EPSJDGetInitialSize(), EPSJDGetKrylovStart()
442: @*/
443: PetscErrorCode EPSJDSetInitialSize(EPS eps,PetscInt initialsize)
444: {

450:   PetscTryMethod(eps,"EPSJDSetInitialSize_C",(EPS,PetscInt),(eps,initialsize));
451:   return(0);
452: }

456: /*@
457:    EPSJDGetFix - Returns the threshold for changing the target in the correction
458:    equation.

460:    Not Collective

462:    Input Parameter:
463: .  eps - the eigenproblem solver context

465:    Output Parameter:
466: .  fix - threshold for changing the target

468:    Note:
469:    The target in the correction equation is fixed at the first iterations.
470:    When the norm of the residual vector is lower than the fix value,
471:    the target is set to the corresponding eigenvalue.

473:    Level: advanced

475: .seealso: EPSJDSetFix()
476: @*/
477: PetscErrorCode EPSJDGetFix(EPS eps,PetscReal *fix)
478: {

484:   PetscUseMethod(eps,"EPSJDGetFix_C",(EPS,PetscReal*),(eps,fix));
485:   return(0);
486: }

490: /*@
491:    EPSJDSetFix - Sets the threshold for changing the target in the correction
492:    equation.

494:    Logically Collective on EPS

496:    Input Parameters:
497: +  eps - the eigenproblem solver context
498: -  fix - threshold for changing the target

500:    Options Database Key:
501: .  -eps_jd_fix - the fix value

503:    Note:
504:    The target in the correction equation is fixed at the first iterations.
505:    When the norm of the residual vector is lower than the fix value,
506:    the target is set to the corresponding eigenvalue.

508:    Level: advanced

510: .seealso: EPSJDGetFix()
511: @*/
512: PetscErrorCode EPSJDSetFix(EPS eps,PetscReal fix)
513: {

519:   PetscTryMethod(eps,"EPSJDSetFix_C",(EPS,PetscReal),(eps,fix));
520:   return(0);
521: }

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

530:    Logically Collective on EPS

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

536:    Options Database Key:
537: .  -eps_jd_const_correction_tol - Deactivates the dynamic stopping criterion.

539:    Level: advanced

541: .seealso: EPSJDGetConstCorrectionTol()
542: @*/
543: PetscErrorCode EPSJDSetConstCorrectionTol(EPS eps,PetscBool constant)
544: {

550:   PetscTryMethod(eps,"EPSJDSetConstCorrectionTol_C",(EPS,PetscBool),(eps,constant));
551:   return(0);
552: }

556: /*@
557:    EPSJDGetConstCorrectionTol - Returns a flag indicating if the dynamic stopping is being used for
558:    solving the correction equation. If the flag is false the KSP relative tolerance is set
559:    to 0.5**i, where i is the number of EPS iterations from the last converged value.

561:    Not Collective

563:    Input Parameter:
564: .  eps - the eigenproblem solver context

566:    Output Parameters:
567: .  constant - boolean flag indicating if the dynamic stopping criterion is not being used.

569:    Level: advanced

571: .seealso: EPSJDGetConstCorrectionTol()
572: @*/
573: PetscErrorCode EPSJDGetConstCorrectionTol(EPS eps,PetscBool *constant)
574: {

580:   PetscUseMethod(eps,"EPSJDGetConstCorrectionTol_C",(EPS,PetscBool*),(eps,constant));
581:   return(0);
582: }

586: /*@
587:    EPSJDGetWindowSizes - Gets the number of converged vectors in the projected
588:    problem (or Rayleigh quotient) and in the projector employed in the correction
589:    equation.

591:    Not Collective

593:    Input Parameter:
594: .  eps - the eigenproblem solver context

596:    Output Parameter:
597: +  pwindow - number of converged vectors in the projector
598: -  qwindow - number of converged vectors in the projected problem

600:    Level: advanced

602: .seealso: EPSJDSetWindowSizes()
603: @*/
604: PetscErrorCode EPSJDGetWindowSizes(EPS eps,PetscInt *pwindow,PetscInt *qwindow)
605: {

610:   PetscUseMethod(eps,"EPSJDGetWindowSizes_C",(EPS,PetscInt*,PetscInt*),(eps,pwindow,qwindow));
611:   return(0);
612: }

616: /*@
617:    EPSJDSetWindowSizes - Sets the number of converged vectors in the projected
618:    problem (or Rayleigh quotient) and in the projector employed in the correction
619:    equation.

621:    Logically Collective on EPS

623:    Input Parameters:
624: +  eps - the eigenproblem solver context
625: .  pwindow - number of converged vectors in the projector
626: -  qwindow - number of converged vectors in the projected problem

628:    Options Database Keys:
629: +  -eps_jd_pwindow - set the number of converged vectors in the projector
630: -  -eps_jd_qwindow - set the number of converged vectors in the projected problem

632:    Level: advanced

634: .seealso: EPSJDGetWindowSizes()
635: @*/
636: PetscErrorCode EPSJDSetWindowSizes(EPS eps,PetscInt pwindow,PetscInt qwindow)
637: {

644:   PetscTryMethod(eps,"EPSJDSetWindowSizes_C",(EPS,PetscInt,PetscInt),(eps,pwindow,qwindow));
645:   return(0);
646: }

650: /*@
651:    EPSJDSetBOrth - Selects the orthogonalization that will be used in the search
652:    subspace in case of generalized Hermitian problems.

654:    Logically Collective on EPS

656:    Input Parameters:
657: +  eps   - the eigenproblem solver context
658: -  borth - whether to B-orthogonalize the search subspace

660:    Options Database Key:
661: .  -eps_jd_borth - Set the orthogonalization used in the search subspace

663:    Level: advanced

665: .seealso: EPSJDGetBOrth()
666: @*/
667: PetscErrorCode EPSJDSetBOrth(EPS eps,PetscBool borth)
668: {

674:   PetscTryMethod(eps,"EPSJDSetBOrth_C",(EPS,PetscBool),(eps,borth));
675:   return(0);
676: }

680: /*@
681:    EPSJDGetBOrth - Returns the orthogonalization used in the search
682:    subspace in case of generalized Hermitian problems.

684:    Not Collective

686:    Input Parameter:
687: .  eps - the eigenproblem solver context

689:    Output Parameters:
690: .  borth - whether to B-orthogonalize the search subspace

692:    Level: advanced

694: .seealso: EPSJDSetBOrth()
695: @*/
696: PetscErrorCode EPSJDGetBOrth(EPS eps,PetscBool *borth)
697: {

703:   PetscUseMethod(eps,"EPSJDGetBOrth_C",(EPS,PetscBool*),(eps,borth));
704:   return(0);
705: }

709: PETSC_EXTERN PetscErrorCode EPSCreate_JD(EPS eps)
710: {
712:   EPS_DAVIDSON   *data;

715:   PetscNewLog(eps,&data);
716:   eps->data = (void*)data;

718:   data->blocksize   = 1;
719:   data->initialsize = 6;
720:   data->minv        = 6;
721:   data->plusk       = 0;
722:   data->ipB         = PETSC_TRUE;
723:   data->fix         = 0.01;
724:   data->krylovstart = PETSC_FALSE;
725:   data->dynamic     = PETSC_FALSE;
726:   data->cX_in_proj  = 0;
727:   data->cX_in_impr  = 0;

729:   eps->ops->solve          = EPSSolve_XD;
730:   eps->ops->setup          = EPSSetUp_XD;
731:   eps->ops->reset          = EPSReset_XD;
732:   eps->ops->backtransform  = EPSBackTransform_Default;
733:   eps->ops->computevectors = EPSComputeVectors_XD;
734:   eps->ops->view           = EPSView_JD;
735:   eps->ops->setfromoptions = EPSSetFromOptions_JD;
736:   eps->ops->setup          = EPSSetUp_JD;
737:   eps->ops->destroy        = EPSDestroy_JD;

739:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",EPSXDSetKrylovStart_XD);
740:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",EPSXDGetKrylovStart_XD);
741:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",EPSXDSetBlockSize_XD);
742:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",EPSXDGetBlockSize_XD);
743:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",EPSXDSetRestart_XD);
744:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",EPSXDGetRestart_XD);
745:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",EPSXDSetInitialSize_XD);
746:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",EPSXDGetInitialSize_XD);
747:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",EPSJDSetFix_JD);
748:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",EPSXDGetFix_XD);
749:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",EPSJDSetConstCorrectionTol_JD);
750:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",EPSJDGetConstCorrectionTol_JD);
751:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetWindowSizes_C",EPSXDSetWindowSizes_XD);
752:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetWindowSizes_C",EPSXDGetWindowSizes_XD);
753:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",EPSXDSetBOrth_XD);
754:   PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",EPSXDGetBOrth_XD);
755:   return(0);
756: }