Actual source code: jd.c
slepc-3.21.0 2024-03-30
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: static PetscErrorCode EPSSetFromOptions_JD(EPS eps,PetscOptionItems *PetscOptionsObject)
35: {
36: PetscBool flg,flg2,op,orth;
37: PetscInt opi,opi0;
38: PetscReal opf;
40: PetscFunctionBegin;
41: PetscOptionsHeadBegin(PetscOptionsObject,"EPS Jacobi-Davidson (JD) Options");
43: PetscCall(EPSJDGetKrylovStart(eps,&op));
44: PetscCall(PetscOptionsBool("-eps_jd_krylov_start","Start the search subspace with a Krylov basis","EPSJDSetKrylovStart",op,&op,&flg));
45: if (flg) PetscCall(EPSJDSetKrylovStart(eps,op));
47: PetscCall(EPSJDGetBOrth(eps,&orth));
48: PetscCall(PetscOptionsBool("-eps_jd_borth","Use B-orthogonalization in the search subspace","EPSJDSetBOrth",op,&op,&flg));
49: if (flg) PetscCall(EPSJDSetBOrth(eps,op));
51: PetscCall(EPSJDGetBlockSize(eps,&opi));
52: PetscCall(PetscOptionsInt("-eps_jd_blocksize","Number of vectors to add to the search subspace","EPSJDSetBlockSize",opi,&opi,&flg));
53: if (flg) PetscCall(EPSJDSetBlockSize(eps,opi));
55: PetscCall(EPSJDGetRestart(eps,&opi,&opi0));
56: PetscCall(PetscOptionsInt("-eps_jd_minv","Size of the search subspace after restarting","EPSJDSetRestart",opi,&opi,&flg));
57: PetscCall(PetscOptionsInt("-eps_jd_plusk","Number of eigenvectors saved from the previous iteration when restarting","EPSJDSetRestart",opi0,&opi0,&flg2));
58: if (flg || flg2) PetscCall(EPSJDSetRestart(eps,opi,opi0));
60: PetscCall(EPSJDGetInitialSize(eps,&opi));
61: PetscCall(PetscOptionsInt("-eps_jd_initial_size","Initial size of the search subspace","EPSJDSetInitialSize",opi,&opi,&flg));
62: if (flg) PetscCall(EPSJDSetInitialSize(eps,opi));
64: PetscCall(EPSJDGetFix(eps,&opf));
65: PetscCall(PetscOptionsReal("-eps_jd_fix","Tolerance for changing the target in the correction equation","EPSJDSetFix",opf,&opf,&flg));
66: if (flg) PetscCall(EPSJDSetFix(eps,opf));
68: PetscCall(EPSJDGetConstCorrectionTol(eps,&op));
69: PetscCall(PetscOptionsBool("-eps_jd_const_correction_tol","Disable the dynamic stopping criterion when solving the correction equation","EPSJDSetConstCorrectionTol",op,&op,&flg));
70: if (flg) PetscCall(EPSJDSetConstCorrectionTol(eps,op));
72: PetscOptionsHeadEnd();
73: PetscFunctionReturn(PETSC_SUCCESS);
74: }
76: static PetscErrorCode EPSSetDefaultST_JD(EPS eps)
77: {
78: KSP ksp;
80: PetscFunctionBegin;
81: if (!((PetscObject)eps->st)->type_name) {
82: PetscCall(STSetType(eps->st,STPRECOND));
83: PetscCall(STPrecondSetKSPHasMat(eps->st,PETSC_TRUE));
84: }
85: PetscCall(STGetKSP(eps->st,&ksp));
86: if (!((PetscObject)ksp)->type_name) {
87: PetscCall(KSPSetType(ksp,KSPBCGSL));
88: PetscCall(KSPSetTolerances(ksp,1e-4,PETSC_DEFAULT,PETSC_DEFAULT,90));
89: }
90: PetscFunctionReturn(PETSC_SUCCESS);
91: }
93: static PetscErrorCode EPSSetUp_JD(EPS eps)
94: {
95: PetscBool t;
96: KSP ksp;
98: PetscFunctionBegin;
99: /* Setup common for all davidson solvers */
100: PetscCall(EPSSetUp_XD(eps));
102: /* Check some constraints */
103: PetscCall(STGetKSP(eps->st,&ksp));
104: PetscCall(PetscObjectTypeCompare((PetscObject)ksp,KSPPREONLY,&t));
105: PetscCheck(!t,PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"EPSJD does not work with KSPPREONLY");
106: PetscFunctionReturn(PETSC_SUCCESS);
107: }
109: static PetscErrorCode EPSView_JD(EPS eps,PetscViewer viewer)
110: {
111: PetscBool isascii,opb;
112: PetscReal opf;
113: PetscInt opi,opi0;
115: PetscFunctionBegin;
116: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
117: if (isascii) {
118: PetscCall(EPSXDGetBOrth_XD(eps,&opb));
119: if (opb) PetscCall(PetscViewerASCIIPrintf(viewer," search subspace is B-orthogonalized\n"));
120: else PetscCall(PetscViewerASCIIPrintf(viewer," search subspace is orthogonalized\n"));
121: PetscCall(EPSXDGetBlockSize_XD(eps,&opi));
122: PetscCall(PetscViewerASCIIPrintf(viewer," block size=%" PetscInt_FMT "\n",opi));
123: PetscCall(EPSXDGetKrylovStart_XD(eps,&opb));
124: if (!opb) PetscCall(PetscViewerASCIIPrintf(viewer," type of the initial subspace: non-Krylov\n"));
125: else PetscCall(PetscViewerASCIIPrintf(viewer," type of the initial subspace: Krylov\n"));
126: PetscCall(EPSXDGetRestart_XD(eps,&opi,&opi0));
127: PetscCall(PetscViewerASCIIPrintf(viewer," size of the subspace after restarting: %" PetscInt_FMT "\n",opi));
128: PetscCall(PetscViewerASCIIPrintf(viewer," number of vectors after restarting from the previous iteration: %" PetscInt_FMT "\n",opi0));
130: PetscCall(EPSJDGetFix_JD(eps,&opf));
131: PetscCall(PetscViewerASCIIPrintf(viewer," threshold for changing the target in the correction equation (fix): %g\n",(double)opf));
133: PetscCall(EPSJDGetConstCorrectionTol_JD(eps,&opb));
134: if (!opb) PetscCall(PetscViewerASCIIPrintf(viewer," using dynamic tolerance for the correction equation\n"));
135: }
136: PetscFunctionReturn(PETSC_SUCCESS);
137: }
139: static PetscErrorCode EPSDestroy_JD(EPS eps)
140: {
141: PetscFunctionBegin;
142: PetscCall(PetscFree(eps->data));
143: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",NULL));
144: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",NULL));
145: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",NULL));
146: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",NULL));
147: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",NULL));
148: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",NULL));
149: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",NULL));
150: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",NULL));
151: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",NULL));
152: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",NULL));
153: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",NULL));
154: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",NULL));
155: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",NULL));
156: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",NULL));
157: PetscFunctionReturn(PETSC_SUCCESS);
158: }
160: /*@
161: EPSJDSetKrylovStart - Activates or deactivates starting the searching
162: subspace with a Krylov basis.
164: Logically Collective
166: Input Parameters:
167: + eps - the eigenproblem solver context
168: - krylovstart - boolean flag
170: Options Database Key:
171: . -eps_jd_krylov_start - Activates starting the searching subspace with a
172: Krylov basis
174: Level: advanced
176: .seealso: EPSJDGetKrylovStart()
177: @*/
178: PetscErrorCode EPSJDSetKrylovStart(EPS eps,PetscBool krylovstart)
179: {
180: PetscFunctionBegin;
183: PetscTryMethod(eps,"EPSJDSetKrylovStart_C",(EPS,PetscBool),(eps,krylovstart));
184: PetscFunctionReturn(PETSC_SUCCESS);
185: }
187: /*@
188: EPSJDGetKrylovStart - Returns a flag indicating if the searching subspace is started with a
189: Krylov basis.
191: Not Collective
193: Input Parameter:
194: . eps - the eigenproblem solver context
196: Output Parameters:
197: . krylovstart - boolean flag indicating if the searching subspace is started
198: with a Krylov basis
200: Level: advanced
202: .seealso: EPSJDSetKrylovStart()
203: @*/
204: PetscErrorCode EPSJDGetKrylovStart(EPS eps,PetscBool *krylovstart)
205: {
206: PetscFunctionBegin;
208: PetscAssertPointer(krylovstart,2);
209: PetscUseMethod(eps,"EPSJDGetKrylovStart_C",(EPS,PetscBool*),(eps,krylovstart));
210: PetscFunctionReturn(PETSC_SUCCESS);
211: }
213: /*@
214: EPSJDSetBlockSize - Sets the number of vectors to be added to the searching space
215: in every iteration.
217: Logically Collective
219: Input Parameters:
220: + eps - the eigenproblem solver context
221: - blocksize - number of vectors added to the search space in every iteration
223: Options Database Key:
224: . -eps_jd_blocksize - number of vectors added to the searching space every iteration
226: Level: advanced
228: .seealso: EPSJDSetKrylovStart()
229: @*/
230: PetscErrorCode EPSJDSetBlockSize(EPS eps,PetscInt blocksize)
231: {
232: PetscFunctionBegin;
235: PetscTryMethod(eps,"EPSJDSetBlockSize_C",(EPS,PetscInt),(eps,blocksize));
236: PetscFunctionReturn(PETSC_SUCCESS);
237: }
239: /*@
240: EPSJDGetBlockSize - Returns the number of vectors to be added to the searching space
241: in every iteration.
243: Not Collective
245: Input Parameter:
246: . eps - the eigenproblem solver context
248: Output Parameter:
249: . blocksize - number of vectors added to the search space in every iteration
251: Level: advanced
253: .seealso: EPSJDSetBlockSize()
254: @*/
255: PetscErrorCode EPSJDGetBlockSize(EPS eps,PetscInt *blocksize)
256: {
257: PetscFunctionBegin;
259: PetscAssertPointer(blocksize,2);
260: PetscUseMethod(eps,"EPSJDGetBlockSize_C",(EPS,PetscInt*),(eps,blocksize));
261: PetscFunctionReturn(PETSC_SUCCESS);
262: }
264: /*@
265: EPSJDSetRestart - Sets the number of vectors of the searching space after
266: restarting and the number of vectors saved from the previous iteration.
268: Logically Collective
270: Input Parameters:
271: + eps - the eigenproblem solver context
272: . minv - number of vectors of the searching subspace after restarting
273: - plusk - number of vectors saved from the previous iteration
275: Options Database Keys:
276: + -eps_jd_minv - number of vectors of the searching subspace after restarting
277: - -eps_jd_plusk - number of vectors saved from the previous iteration
279: Level: advanced
281: .seealso: EPSJDGetRestart()
282: @*/
283: PetscErrorCode EPSJDSetRestart(EPS eps,PetscInt minv,PetscInt plusk)
284: {
285: PetscFunctionBegin;
289: PetscTryMethod(eps,"EPSJDSetRestart_C",(EPS,PetscInt,PetscInt),(eps,minv,plusk));
290: PetscFunctionReturn(PETSC_SUCCESS);
291: }
293: /*@
294: EPSJDGetRestart - Gets the number of vectors of the searching space after
295: restarting and the number of vectors saved from the previous iteration.
297: Not Collective
299: Input Parameter:
300: . eps - the eigenproblem solver context
302: Output Parameters:
303: + minv - number of vectors of the searching subspace after restarting
304: - plusk - number of vectors saved from the previous iteration
306: Level: advanced
308: .seealso: EPSJDSetRestart()
309: @*/
310: PetscErrorCode EPSJDGetRestart(EPS eps,PetscInt *minv,PetscInt *plusk)
311: {
312: PetscFunctionBegin;
314: PetscUseMethod(eps,"EPSJDGetRestart_C",(EPS,PetscInt*,PetscInt*),(eps,minv,plusk));
315: PetscFunctionReturn(PETSC_SUCCESS);
316: }
318: /*@
319: EPSJDSetInitialSize - Sets the initial size of the searching space.
321: Logically Collective
323: Input Parameters:
324: + eps - the eigenproblem solver context
325: - initialsize - number of vectors of the initial searching subspace
327: Options Database Key:
328: . -eps_jd_initial_size - number of vectors of the initial searching subspace
330: Notes:
331: If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
332: EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
333: provided vectors are not enough, the solver completes the subspace with
334: random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
335: gets the first vector provided by the user or, if not available, a random vector,
336: and expands the Krylov basis up to initialsize vectors.
338: Level: advanced
340: .seealso: EPSJDGetInitialSize(), EPSJDGetKrylovStart()
341: @*/
342: PetscErrorCode EPSJDSetInitialSize(EPS eps,PetscInt initialsize)
343: {
344: PetscFunctionBegin;
347: PetscTryMethod(eps,"EPSJDSetInitialSize_C",(EPS,PetscInt),(eps,initialsize));
348: PetscFunctionReturn(PETSC_SUCCESS);
349: }
351: /*@
352: EPSJDGetInitialSize - Returns the initial size of the searching space.
354: Not Collective
356: Input Parameter:
357: . eps - the eigenproblem solver context
359: Output Parameter:
360: . initialsize - number of vectors of the initial searching subspace
362: Notes:
363: If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
364: EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
365: provided vectors are not enough, the solver completes the subspace with
366: random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
367: gets the first vector provided by the user or, if not available, a random vector,
368: and expands the Krylov basis up to initialsize vectors.
370: Level: advanced
372: .seealso: EPSJDSetInitialSize(), EPSJDGetKrylovStart()
373: @*/
374: PetscErrorCode EPSJDGetInitialSize(EPS eps,PetscInt *initialsize)
375: {
376: PetscFunctionBegin;
378: PetscAssertPointer(initialsize,2);
379: PetscUseMethod(eps,"EPSJDGetInitialSize_C",(EPS,PetscInt*),(eps,initialsize));
380: PetscFunctionReturn(PETSC_SUCCESS);
381: }
383: static PetscErrorCode EPSJDSetFix_JD(EPS eps,PetscReal fix)
384: {
385: EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;
387: PetscFunctionBegin;
388: if (fix == (PetscReal)PETSC_DEFAULT || fix == (PetscReal)PETSC_DECIDE) fix = 0.01;
389: PetscCheck(fix>0.0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid fix value, must be >0");
390: data->fix = fix;
391: PetscFunctionReturn(PETSC_SUCCESS);
392: }
394: /*@
395: EPSJDSetFix - Sets the threshold for changing the target in the correction
396: equation.
398: Logically Collective
400: Input Parameters:
401: + eps - the eigenproblem solver context
402: - fix - threshold for changing the target
404: Options Database Key:
405: . -eps_jd_fix - the fix value
407: Note:
408: The target in the correction equation is fixed at the first iterations.
409: When the norm of the residual vector is lower than the fix value,
410: the target is set to the corresponding eigenvalue.
412: Level: advanced
414: .seealso: EPSJDGetFix()
415: @*/
416: PetscErrorCode EPSJDSetFix(EPS eps,PetscReal fix)
417: {
418: PetscFunctionBegin;
421: PetscTryMethod(eps,"EPSJDSetFix_C",(EPS,PetscReal),(eps,fix));
422: PetscFunctionReturn(PETSC_SUCCESS);
423: }
425: PetscErrorCode EPSJDGetFix_JD(EPS eps,PetscReal *fix)
426: {
427: EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;
429: PetscFunctionBegin;
430: *fix = data->fix;
431: PetscFunctionReturn(PETSC_SUCCESS);
432: }
434: /*@
435: EPSJDGetFix - Returns the threshold for changing the target in the correction
436: equation.
438: Not Collective
440: Input Parameter:
441: . eps - the eigenproblem solver context
443: Output Parameter:
444: . fix - threshold for changing the target
446: Note:
447: The target in the correction equation is fixed at the first iterations.
448: When the norm of the residual vector is lower than the fix value,
449: the target is set to the corresponding eigenvalue.
451: Level: advanced
453: .seealso: EPSJDSetFix()
454: @*/
455: PetscErrorCode EPSJDGetFix(EPS eps,PetscReal *fix)
456: {
457: PetscFunctionBegin;
459: PetscAssertPointer(fix,2);
460: PetscUseMethod(eps,"EPSJDGetFix_C",(EPS,PetscReal*),(eps,fix));
461: PetscFunctionReturn(PETSC_SUCCESS);
462: }
464: static PetscErrorCode EPSJDSetConstCorrectionTol_JD(EPS eps,PetscBool constant)
465: {
466: EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;
468: PetscFunctionBegin;
469: data->dynamic = PetscNot(constant);
470: PetscFunctionReturn(PETSC_SUCCESS);
471: }
473: /*@
474: EPSJDSetConstCorrectionTol - If true, deactivates the dynamic stopping criterion
475: (also called Newton) that sets the KSP relative tolerance
476: to 0.5**i, where i is the number of EPS iterations from the last converged value.
478: Logically Collective
480: Input Parameters:
481: + eps - the eigenproblem solver context
482: - constant - if false, the KSP relative tolerance is set to 0.5**i.
484: Options Database Key:
485: . -eps_jd_const_correction_tol - Deactivates the dynamic stopping criterion.
487: Level: advanced
489: .seealso: EPSJDGetConstCorrectionTol()
490: @*/
491: PetscErrorCode EPSJDSetConstCorrectionTol(EPS eps,PetscBool constant)
492: {
493: PetscFunctionBegin;
496: PetscTryMethod(eps,"EPSJDSetConstCorrectionTol_C",(EPS,PetscBool),(eps,constant));
497: PetscFunctionReturn(PETSC_SUCCESS);
498: }
500: PetscErrorCode EPSJDGetConstCorrectionTol_JD(EPS eps,PetscBool *constant)
501: {
502: EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;
504: PetscFunctionBegin;
505: *constant = PetscNot(data->dynamic);
506: PetscFunctionReturn(PETSC_SUCCESS);
507: }
509: /*@
510: EPSJDGetConstCorrectionTol - Returns a flag indicating if the dynamic stopping is being used for
511: solving the correction equation.
513: Not Collective
515: Input Parameter:
516: . eps - the eigenproblem solver context
518: Output Parameter:
519: . constant - boolean flag indicating if the dynamic stopping criterion is not being used.
521: Notes:
522: If the flag is false the KSP relative tolerance is set to 0.5**i, where i is the number
523: of EPS iterations from the last converged value.
525: Level: advanced
527: .seealso: EPSJDSetConstCorrectionTol()
528: @*/
529: PetscErrorCode EPSJDGetConstCorrectionTol(EPS eps,PetscBool *constant)
530: {
531: PetscFunctionBegin;
533: PetscAssertPointer(constant,2);
534: PetscUseMethod(eps,"EPSJDGetConstCorrectionTol_C",(EPS,PetscBool*),(eps,constant));
535: PetscFunctionReturn(PETSC_SUCCESS);
536: }
538: /*@
539: EPSJDSetBOrth - Selects the orthogonalization that will be used in the search
540: subspace in case of generalized Hermitian problems.
542: Logically Collective
544: Input Parameters:
545: + eps - the eigenproblem solver context
546: - borth - whether to B-orthogonalize the search subspace
548: Options Database Key:
549: . -eps_jd_borth - Set the orthogonalization used in the search subspace
551: Level: advanced
553: .seealso: EPSJDGetBOrth()
554: @*/
555: PetscErrorCode EPSJDSetBOrth(EPS eps,PetscBool borth)
556: {
557: PetscFunctionBegin;
560: PetscTryMethod(eps,"EPSJDSetBOrth_C",(EPS,PetscBool),(eps,borth));
561: PetscFunctionReturn(PETSC_SUCCESS);
562: }
564: /*@
565: EPSJDGetBOrth - Returns the orthogonalization used in the search
566: subspace in case of generalized Hermitian problems.
568: Not Collective
570: Input Parameter:
571: . eps - the eigenproblem solver context
573: Output Parameters:
574: . borth - whether to B-orthogonalize the search subspace
576: Level: advanced
578: .seealso: EPSJDSetBOrth()
579: @*/
580: PetscErrorCode EPSJDGetBOrth(EPS eps,PetscBool *borth)
581: {
582: PetscFunctionBegin;
584: PetscAssertPointer(borth,2);
585: PetscUseMethod(eps,"EPSJDGetBOrth_C",(EPS,PetscBool*),(eps,borth));
586: PetscFunctionReturn(PETSC_SUCCESS);
587: }
589: SLEPC_EXTERN PetscErrorCode EPSCreate_JD(EPS eps)
590: {
591: EPS_DAVIDSON *data;
593: PetscFunctionBegin;
594: PetscCall(PetscNew(&data));
595: eps->data = (void*)data;
597: data->blocksize = 1;
598: data->initialsize = 0;
599: data->minv = 0;
600: data->plusk = PETSC_DEFAULT;
601: data->ipB = PETSC_TRUE;
602: data->fix = 0.01;
603: data->krylovstart = PETSC_FALSE;
604: data->dynamic = PETSC_FALSE;
606: eps->useds = PETSC_TRUE;
607: eps->categ = EPS_CATEGORY_PRECOND;
609: eps->ops->solve = EPSSolve_XD;
610: eps->ops->setup = EPSSetUp_JD;
611: eps->ops->setupsort = EPSSetUpSort_Default;
612: eps->ops->setfromoptions = EPSSetFromOptions_JD;
613: eps->ops->destroy = EPSDestroy_JD;
614: eps->ops->reset = EPSReset_XD;
615: eps->ops->view = EPSView_JD;
616: eps->ops->backtransform = EPSBackTransform_Default;
617: eps->ops->computevectors = EPSComputeVectors_XD;
618: eps->ops->setdefaultst = EPSSetDefaultST_JD;
620: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",EPSXDSetKrylovStart_XD));
621: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",EPSXDGetKrylovStart_XD));
622: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",EPSXDSetBlockSize_XD));
623: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",EPSXDGetBlockSize_XD));
624: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",EPSXDSetRestart_XD));
625: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",EPSXDGetRestart_XD));
626: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",EPSXDSetInitialSize_XD));
627: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",EPSXDGetInitialSize_XD));
628: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",EPSJDSetFix_JD));
629: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",EPSJDGetFix_JD));
630: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",EPSJDSetConstCorrectionTol_JD));
631: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",EPSJDGetConstCorrectionTol_JD));
632: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",EPSXDSetBOrth_XD));
633: PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",EPSXDGetBOrth_XD));
634: PetscFunctionReturn(PETSC_SUCCESS);
635: }