Actual source code: jd.c

slepc-main 2024-11-09
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: 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_CURRENT,PETSC_CURRENT,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:    Note:
280:    PETSC_CURRENT can be used to preserve the current value of any of the
281:    arguments, and PETSC_DETERMINE to set them to a default value.

283:    Level: advanced

285: .seealso: EPSJDGetRestart()
286: @*/
287: PetscErrorCode EPSJDSetRestart(EPS eps,PetscInt minv,PetscInt plusk)
288: {
289:   PetscFunctionBegin;
293:   PetscTryMethod(eps,"EPSJDSetRestart_C",(EPS,PetscInt,PetscInt),(eps,minv,plusk));
294:   PetscFunctionReturn(PETSC_SUCCESS);
295: }

297: /*@
298:    EPSJDGetRestart - Gets the number of vectors of the searching space after
299:    restarting and the number of vectors saved from the previous iteration.

301:    Not Collective

303:    Input Parameter:
304: .  eps - the eigenproblem solver context

306:    Output Parameters:
307: +  minv - number of vectors of the searching subspace after restarting
308: -  plusk - number of vectors saved from the previous iteration

310:    Level: advanced

312: .seealso: EPSJDSetRestart()
313: @*/
314: PetscErrorCode EPSJDGetRestart(EPS eps,PetscInt *minv,PetscInt *plusk)
315: {
316:   PetscFunctionBegin;
318:   PetscUseMethod(eps,"EPSJDGetRestart_C",(EPS,PetscInt*,PetscInt*),(eps,minv,plusk));
319:   PetscFunctionReturn(PETSC_SUCCESS);
320: }

322: /*@
323:    EPSJDSetInitialSize - Sets the initial size of the searching space.

325:    Logically Collective

327:    Input Parameters:
328: +  eps - the eigenproblem solver context
329: -  initialsize - number of vectors of the initial searching subspace

331:    Options Database Key:
332: .  -eps_jd_initial_size - number of vectors of the initial searching subspace

334:    Notes:
335:    If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
336:    EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
337:    provided vectors are not enough, the solver completes the subspace with
338:    random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
339:    gets the first vector provided by the user or, if not available, a random vector,
340:    and expands the Krylov basis up to initialsize vectors.

342:    Level: advanced

344: .seealso: EPSJDGetInitialSize(), EPSJDGetKrylovStart()
345: @*/
346: PetscErrorCode EPSJDSetInitialSize(EPS eps,PetscInt initialsize)
347: {
348:   PetscFunctionBegin;
351:   PetscTryMethod(eps,"EPSJDSetInitialSize_C",(EPS,PetscInt),(eps,initialsize));
352:   PetscFunctionReturn(PETSC_SUCCESS);
353: }

355: /*@
356:    EPSJDGetInitialSize - Returns the initial size of the searching space.

358:    Not Collective

360:    Input Parameter:
361: .  eps - the eigenproblem solver context

363:    Output Parameter:
364: .  initialsize - number of vectors of the initial searching subspace

366:    Notes:
367:    If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
368:    EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
369:    provided vectors are not enough, the solver completes the subspace with
370:    random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
371:    gets the first vector provided by the user or, if not available, a random vector,
372:    and expands the Krylov basis up to initialsize vectors.

374:    Level: advanced

376: .seealso: EPSJDSetInitialSize(), EPSJDGetKrylovStart()
377: @*/
378: PetscErrorCode EPSJDGetInitialSize(EPS eps,PetscInt *initialsize)
379: {
380:   PetscFunctionBegin;
382:   PetscAssertPointer(initialsize,2);
383:   PetscUseMethod(eps,"EPSJDGetInitialSize_C",(EPS,PetscInt*),(eps,initialsize));
384:   PetscFunctionReturn(PETSC_SUCCESS);
385: }

387: static PetscErrorCode EPSJDSetFix_JD(EPS eps,PetscReal fix)
388: {
389:   EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;

391:   PetscFunctionBegin;
392:   if (fix == (PetscReal)PETSC_DEFAULT || fix == (PetscReal)PETSC_DECIDE) fix = 0.01;
393:   PetscCheck(fix>0.0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid fix value, must be >0");
394:   data->fix = fix;
395:   PetscFunctionReturn(PETSC_SUCCESS);
396: }

398: /*@
399:    EPSJDSetFix - Sets the threshold for changing the target in the correction
400:    equation.

402:    Logically Collective

404:    Input Parameters:
405: +  eps - the eigenproblem solver context
406: -  fix - threshold for changing the target

408:    Options Database Key:
409: .  -eps_jd_fix - the fix value

411:    Note:
412:    The target in the correction equation is fixed at the first iterations.
413:    When the norm of the residual vector is lower than the fix value,
414:    the target is set to the corresponding eigenvalue.

416:    Level: advanced

418: .seealso: EPSJDGetFix()
419: @*/
420: PetscErrorCode EPSJDSetFix(EPS eps,PetscReal fix)
421: {
422:   PetscFunctionBegin;
425:   PetscTryMethod(eps,"EPSJDSetFix_C",(EPS,PetscReal),(eps,fix));
426:   PetscFunctionReturn(PETSC_SUCCESS);
427: }

429: PetscErrorCode EPSJDGetFix_JD(EPS eps,PetscReal *fix)
430: {
431:   EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;

433:   PetscFunctionBegin;
434:   *fix = data->fix;
435:   PetscFunctionReturn(PETSC_SUCCESS);
436: }

438: /*@
439:    EPSJDGetFix - Returns the threshold for changing the target in the correction
440:    equation.

442:    Not Collective

444:    Input Parameter:
445: .  eps - the eigenproblem solver context

447:    Output Parameter:
448: .  fix - threshold for changing the target

450:    Note:
451:    The target in the correction equation is fixed at the first iterations.
452:    When the norm of the residual vector is lower than the fix value,
453:    the target is set to the corresponding eigenvalue.

455:    Level: advanced

457: .seealso: EPSJDSetFix()
458: @*/
459: PetscErrorCode EPSJDGetFix(EPS eps,PetscReal *fix)
460: {
461:   PetscFunctionBegin;
463:   PetscAssertPointer(fix,2);
464:   PetscUseMethod(eps,"EPSJDGetFix_C",(EPS,PetscReal*),(eps,fix));
465:   PetscFunctionReturn(PETSC_SUCCESS);
466: }

468: static PetscErrorCode EPSJDSetConstCorrectionTol_JD(EPS eps,PetscBool constant)
469: {
470:   EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;

472:   PetscFunctionBegin;
473:   data->dynamic = PetscNot(constant);
474:   PetscFunctionReturn(PETSC_SUCCESS);
475: }

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

482:    Logically Collective

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

488:    Options Database Key:
489: .  -eps_jd_const_correction_tol - Deactivates the dynamic stopping criterion.

491:    Level: advanced

493: .seealso: EPSJDGetConstCorrectionTol()
494: @*/
495: PetscErrorCode EPSJDSetConstCorrectionTol(EPS eps,PetscBool constant)
496: {
497:   PetscFunctionBegin;
500:   PetscTryMethod(eps,"EPSJDSetConstCorrectionTol_C",(EPS,PetscBool),(eps,constant));
501:   PetscFunctionReturn(PETSC_SUCCESS);
502: }

504: PetscErrorCode EPSJDGetConstCorrectionTol_JD(EPS eps,PetscBool *constant)
505: {
506:   EPS_DAVIDSON *data = (EPS_DAVIDSON*)eps->data;

508:   PetscFunctionBegin;
509:   *constant = PetscNot(data->dynamic);
510:   PetscFunctionReturn(PETSC_SUCCESS);
511: }

513: /*@
514:    EPSJDGetConstCorrectionTol - Returns a flag indicating if the dynamic stopping is being used for
515:    solving the correction equation.

517:    Not Collective

519:    Input Parameter:
520: .  eps - the eigenproblem solver context

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

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

529:    Level: advanced

531: .seealso: EPSJDSetConstCorrectionTol()
532: @*/
533: PetscErrorCode EPSJDGetConstCorrectionTol(EPS eps,PetscBool *constant)
534: {
535:   PetscFunctionBegin;
537:   PetscAssertPointer(constant,2);
538:   PetscUseMethod(eps,"EPSJDGetConstCorrectionTol_C",(EPS,PetscBool*),(eps,constant));
539:   PetscFunctionReturn(PETSC_SUCCESS);
540: }

542: /*@
543:    EPSJDSetBOrth - Selects the orthogonalization that will be used in the search
544:    subspace in case of generalized Hermitian problems.

546:    Logically Collective

548:    Input Parameters:
549: +  eps   - the eigenproblem solver context
550: -  borth - whether to B-orthogonalize the search subspace

552:    Options Database Key:
553: .  -eps_jd_borth - Set the orthogonalization used in the search subspace

555:    Level: advanced

557: .seealso: EPSJDGetBOrth()
558: @*/
559: PetscErrorCode EPSJDSetBOrth(EPS eps,PetscBool borth)
560: {
561:   PetscFunctionBegin;
564:   PetscTryMethod(eps,"EPSJDSetBOrth_C",(EPS,PetscBool),(eps,borth));
565:   PetscFunctionReturn(PETSC_SUCCESS);
566: }

568: /*@
569:    EPSJDGetBOrth - Returns the orthogonalization used in the search
570:    subspace in case of generalized Hermitian problems.

572:    Not Collective

574:    Input Parameter:
575: .  eps - the eigenproblem solver context

577:    Output Parameters:
578: .  borth - whether to B-orthogonalize the search subspace

580:    Level: advanced

582: .seealso: EPSJDSetBOrth()
583: @*/
584: PetscErrorCode EPSJDGetBOrth(EPS eps,PetscBool *borth)
585: {
586:   PetscFunctionBegin;
588:   PetscAssertPointer(borth,2);
589:   PetscUseMethod(eps,"EPSJDGetBOrth_C",(EPS,PetscBool*),(eps,borth));
590:   PetscFunctionReturn(PETSC_SUCCESS);
591: }

593: SLEPC_EXTERN PetscErrorCode EPSCreate_JD(EPS eps)
594: {
595:   EPS_DAVIDSON   *data;

597:   PetscFunctionBegin;
598:   PetscCall(PetscNew(&data));
599:   eps->data = (void*)data;

601:   data->blocksize   = 1;
602:   data->initialsize = 0;
603:   data->minv        = 0;
604:   data->plusk       = PETSC_DETERMINE;
605:   data->ipB         = PETSC_TRUE;
606:   data->fix         = 0.01;
607:   data->krylovstart = PETSC_FALSE;
608:   data->dynamic     = PETSC_FALSE;

610:   eps->useds = PETSC_TRUE;
611:   eps->categ = EPS_CATEGORY_PRECOND;

613:   eps->ops->solve          = EPSSolve_XD;
614:   eps->ops->setup          = EPSSetUp_JD;
615:   eps->ops->setupsort      = EPSSetUpSort_Default;
616:   eps->ops->setfromoptions = EPSSetFromOptions_JD;
617:   eps->ops->destroy        = EPSDestroy_JD;
618:   eps->ops->reset          = EPSReset_XD;
619:   eps->ops->view           = EPSView_JD;
620:   eps->ops->backtransform  = EPSBackTransform_Default;
621:   eps->ops->computevectors = EPSComputeVectors_XD;
622:   eps->ops->setdefaultst   = EPSSetDefaultST_JD;

624:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",EPSXDSetKrylovStart_XD));
625:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",EPSXDGetKrylovStart_XD));
626:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",EPSXDSetBlockSize_XD));
627:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",EPSXDGetBlockSize_XD));
628:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",EPSXDSetRestart_XD));
629:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",EPSXDGetRestart_XD));
630:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",EPSXDSetInitialSize_XD));
631:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",EPSXDGetInitialSize_XD));
632:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",EPSJDSetFix_JD));
633:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",EPSJDGetFix_JD));
634:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",EPSJDSetConstCorrectionTol_JD));
635:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",EPSJDGetConstCorrectionTol_JD));
636:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",EPSXDSetBOrth_XD));
637:   PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",EPSXDGetBOrth_XD));
638:   PetscFunctionReturn(PETSC_SUCCESS);
639: }