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-2015, 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(PetscOptions *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 EPSDestroy_JD(EPS eps)142: {
143: PetscErrorCode ierr;
146: PetscFree(eps->data);
147: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",NULL);
148: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",NULL);
149: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",NULL);
150: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",NULL);
151: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",NULL);
152: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",NULL);
153: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",NULL);
154: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",NULL);
155: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",NULL);
156: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",NULL);
157: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",NULL);
158: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",NULL);
159: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetWindowSizes_C",NULL);
160: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetWindowSizes_C",NULL);
161: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",NULL);
162: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",NULL);
163: return(0);
164: }
168: /*@
169: EPSJDSetKrylovStart - Activates or deactivates starting the searching
170: subspace with a Krylov basis.
172: Logically Collective on EPS174: Input Parameters:
175: + eps - the eigenproblem solver context
176: - krylovstart - boolean flag
178: Options Database Key:
179: . -eps_jd_krylov_start - Activates starting the searching subspace with a
180: Krylov basis
182: Level: advanced
184: .seealso: EPSJDGetKrylovStart()
185: @*/
186: PetscErrorCode EPSJDSetKrylovStart(EPS eps,PetscBool krylovstart)187: {
193: PetscTryMethod(eps,"EPSJDSetKrylovStart_C",(EPS,PetscBool),(eps,krylovstart));
194: return(0);
195: }
199: /*@
200: EPSJDGetKrylovStart - Returns a flag indicating if the searching subspace is started with a
201: Krylov basis.
203: Not Collective
205: Input Parameter:
206: . eps - the eigenproblem solver context
208: Output Parameters:
209: . krylovstart - boolean flag indicating if the searching subspace is started
210: with a Krylov basis
212: Level: advanced
214: .seealso: EPSJDGetKrylovStart()
215: @*/
216: PetscErrorCode EPSJDGetKrylovStart(EPS eps,PetscBool *krylovstart)217: {
223: PetscTryMethod(eps,"EPSJDGetKrylovStart_C",(EPS,PetscBool*),(eps,krylovstart));
224: return(0);
225: }
229: /*@
230: EPSJDSetBlockSize - Sets the number of vectors to be added to the searching space
231: in every iteration.
233: Logically Collective on EPS235: Input Parameters:
236: + eps - the eigenproblem solver context
237: - blocksize - number of vectors added to the search space in every iteration
239: Options Database Key:
240: . -eps_jd_blocksize - number of vectors added to the searching space every iteration
242: Level: advanced
244: .seealso: EPSJDSetKrylovStart()
245: @*/
246: PetscErrorCode EPSJDSetBlockSize(EPS eps,PetscInt blocksize)247: {
253: PetscTryMethod(eps,"EPSJDSetBlockSize_C",(EPS,PetscInt),(eps,blocksize));
254: return(0);
255: }
259: /*@
260: EPSJDGetBlockSize - Returns the number of vectors to be added to the searching space
261: in every iteration.
263: Not Collective
265: Input Parameter:
266: . eps - the eigenproblem solver context
268: Output Parameter:
269: . blocksize - number of vectors added to the search space in every iteration
271: Level: advanced
273: .seealso: EPSJDSetBlockSize()
274: @*/
275: PetscErrorCode EPSJDGetBlockSize(EPS eps,PetscInt *blocksize)276: {
282: PetscTryMethod(eps,"EPSJDGetBlockSize_C",(EPS,PetscInt*),(eps,blocksize));
283: return(0);
284: }
288: /*@
289: EPSJDGetRestart - Gets the number of vectors of the searching space after
290: restarting and the number of vectors saved from the previous iteration.
292: Not Collective
294: Input Parameter:
295: . eps - the eigenproblem solver context
297: Output Parameter:
298: + minv - number of vectors of the searching subspace after restarting
299: - plusk - number of vectors saved from the previous iteration
301: Level: advanced
303: .seealso: EPSJDSetRestart()
304: @*/
305: PetscErrorCode EPSJDGetRestart(EPS eps,PetscInt *minv,PetscInt *plusk)306: {
311: PetscTryMethod(eps,"EPSJDGetRestart_C",(EPS,PetscInt*,PetscInt*),(eps,minv,plusk));
312: return(0);
313: }
317: /*@
318: EPSJDSetRestart - Sets the number of vectors of the searching space after
319: restarting and the number of vectors saved from the previous iteration.
321: Logically Collective on EPS323: Input Parameters:
324: + eps - the eigenproblem solver context
325: . minv - number of vectors of the searching subspace after restarting
326: - plusk - number of vectors saved from the previous iteration
328: Options Database Keys:
329: + -eps_jd_minv - number of vectors of the searching subspace after restarting
330: - -eps_jd_plusk - number of vectors saved from the previous iteration
332: Level: advanced
334: .seealso: EPSJDGetRestart()
335: @*/
336: PetscErrorCode EPSJDSetRestart(EPS eps,PetscInt minv,PetscInt plusk)337: {
344: PetscTryMethod(eps,"EPSJDSetRestart_C",(EPS,PetscInt,PetscInt),(eps,minv,plusk));
345: return(0);
346: }
350: /*@
351: EPSJDGetInitialSize - Returns the initial size of the searching space.
353: Not Collective
355: Input Parameter:
356: . eps - the eigenproblem solver context
358: Output Parameter:
359: . initialsize - number of vectors of the initial searching subspace
361: Notes:
362: If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
363: EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
364: provided vectors are not enough, the solver completes the subspace with
365: random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
366: gets the first vector provided by the user or, if not available, a random vector,
367: and expands the Krylov basis up to initialsize vectors.
369: Level: advanced
371: .seealso: EPSJDSetInitialSize(), EPSJDGetKrylovStart()
372: @*/
373: PetscErrorCode EPSJDGetInitialSize(EPS eps,PetscInt *initialsize)374: {
380: PetscTryMethod(eps,"EPSJDGetInitialSize_C",(EPS,PetscInt*),(eps,initialsize));
381: return(0);
382: }
386: /*@
387: EPSJDSetInitialSize - Sets the initial size of the searching space.
389: Logically Collective on EPS391: Input Parameters:
392: + eps - the eigenproblem solver context
393: - initialsize - number of vectors of the initial searching subspace
395: Options Database Key:
396: . -eps_jd_initial_size - number of vectors of the initial searching subspace
398: Notes:
399: If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
400: EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
401: provided vectors are not enough, the solver completes the subspace with
402: random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
403: gets the first vector provided by the user or, if not available, a random vector,
404: and expands the Krylov basis up to initialsize vectors.
406: Level: advanced
408: .seealso: EPSJDGetInitialSize(), EPSJDGetKrylovStart()
409: @*/
410: PetscErrorCode EPSJDSetInitialSize(EPS eps,PetscInt initialsize)411: {
417: PetscTryMethod(eps,"EPSJDSetInitialSize_C",(EPS,PetscInt),(eps,initialsize));
418: return(0);
419: }
423: /*@
424: EPSJDGetFix - Returns the threshold for changing the target in the correction
425: equation.
427: Not Collective
429: Input Parameter:
430: . eps - the eigenproblem solver context
432: Output Parameter:
433: . fix - threshold for changing the target
435: Note:
436: The target in the correction equation is fixed at the first iterations.
437: When the norm of the residual vector is lower than the fix value,
438: the target is set to the corresponding eigenvalue.
440: Level: advanced
442: .seealso: EPSJDSetFix()
443: @*/
444: PetscErrorCode EPSJDGetFix(EPS eps,PetscReal *fix)445: {
451: PetscTryMethod(eps,"EPSJDGetFix_C",(EPS,PetscReal*),(eps,fix));
452: return(0);
453: }
457: /*@
458: EPSJDSetFix - Sets the threshold for changing the target in the correction
459: equation.
461: Logically Collective on EPS463: Input Parameters:
464: + eps - the eigenproblem solver context
465: - fix - threshold for changing the target
467: Options Database Key:
468: . -eps_jd_fix - the fix value
470: Note:
471: The target in the correction equation is fixed at the first iterations.
472: When the norm of the residual vector is lower than the fix value,
473: the target is set to the corresponding eigenvalue.
475: Level: advanced
477: .seealso: EPSJDGetFix()
478: @*/
479: PetscErrorCode EPSJDSetFix(EPS eps,PetscReal fix)480: {
486: PetscTryMethod(eps,"EPSJDSetFix_C",(EPS,PetscReal),(eps,fix));
487: return(0);
488: }
492: /*@
493: EPSJDSetConstCorrectionTol - If true, deactivates the dynamic stopping criterion
494: (also called Newton) that sets the KSP relative tolerance
495: to 0.5**i, where i is the number of EPS iterations from the last converged value.
497: Logically Collective on EPS499: Input Parameters:
500: + eps - the eigenproblem solver context
501: - constant - if false, the KSP relative tolerance is set to 0.5**i.
503: Options Database Key:
504: . -eps_jd_const_correction_tol - Deactivates the dynamic stopping criterion.
506: Level: advanced
508: .seealso: EPSJDGetConstCorrectionTol()
509: @*/
510: PetscErrorCode EPSJDSetConstCorrectionTol(EPS eps,PetscBool constant)511: {
517: PetscTryMethod(eps,"EPSJDSetConstCorrectionTol_C",(EPS,PetscBool),(eps,constant));
518: return(0);
519: }
523: /*@
524: EPSJDGetConstCorrectionTol - Returns a flag indicating if the dynamic stopping is being used for
525: solving the correction equation. If the flag is false the KSP relative tolerance is set
526: to 0.5**i, where i is the number of EPS iterations from the last converged value.
528: Not Collective
530: Input Parameter:
531: . eps - the eigenproblem solver context
533: Output Parameters:
534: . constant - boolean flag indicating if the dynamic stopping criterion is not being used.
536: Level: advanced
538: .seealso: EPSJDGetConstCorrectionTol()
539: @*/
540: PetscErrorCode EPSJDGetConstCorrectionTol(EPS eps,PetscBool *constant)541: {
547: PetscTryMethod(eps,"EPSJDGetConstCorrectionTol",(EPS,PetscBool*),(eps,constant));
548: return(0);
549: }
553: /*@
554: EPSJDGetWindowSizes - Gets the number of converged vectors in the projected
555: problem (or Rayleigh quotient) and in the projector employed in the correction
556: equation.
558: Not Collective
560: Input Parameter:
561: . eps - the eigenproblem solver context
563: Output Parameter:
564: + pwindow - number of converged vectors in the projector
565: - qwindow - number of converged vectors in the projected problem
567: Level: advanced
569: .seealso: EPSJDSetWindowSizes()
570: @*/
571: PetscErrorCode EPSJDGetWindowSizes(EPS eps,PetscInt *pwindow,PetscInt *qwindow)572: {
577: PetscTryMethod(eps,"EPSJDGetWindowSizes_C",(EPS,PetscInt*,PetscInt*),(eps,pwindow,qwindow));
578: return(0);
579: }
583: /*@
584: EPSJDSetWindowSizes - Sets the number of converged vectors in the projected
585: problem (or Rayleigh quotient) and in the projector employed in the correction
586: equation.
588: Logically Collective on EPS590: Input Parameters:
591: + eps - the eigenproblem solver context
592: . pwindow - number of converged vectors in the projector
593: - qwindow - number of converged vectors in the projected problem
595: Options Database Keys:
596: + -eps_jd_pwindow - set the number of converged vectors in the projector
597: - -eps_jd_qwindow - set the number of converged vectors in the projected problem
599: Level: advanced
601: .seealso: EPSJDGetWindowSizes()
602: @*/
603: PetscErrorCode EPSJDSetWindowSizes(EPS eps,PetscInt pwindow,PetscInt qwindow)604: {
611: PetscTryMethod(eps,"EPSJDSetWindowSizes_C",(EPS,PetscInt,PetscInt),(eps,pwindow,qwindow));
612: return(0);
613: }
617: /*@
618: EPSJDSetBOrth - Selects the orthogonalization that will be used in the search
619: subspace in case of generalized Hermitian problems.
621: Logically Collective on EPS623: Input Parameters:
624: + eps - the eigenproblem solver context
625: - borth - whether to B-orthogonalize the search subspace
627: Options Database Key:
628: . -eps_jd_borth - Set the orthogonalization used in the search subspace
630: Level: advanced
632: .seealso: EPSJDGetBOrth()
633: @*/
634: PetscErrorCode EPSJDSetBOrth(EPS eps,PetscBool borth)635: {
641: PetscTryMethod(eps,"EPSJDSetBOrth_C",(EPS,PetscBool),(eps,borth));
642: return(0);
643: }
647: /*@
648: EPSJDGetBOrth - Returns the orthogonalization used in the search
649: subspace in case of generalized Hermitian problems.
651: Not Collective
653: Input Parameter:
654: . eps - the eigenproblem solver context
656: Output Parameters:
657: . borth - whether to B-orthogonalize the search subspace
659: Level: advanced
661: .seealso: EPSJDSetBOrth()
662: @*/
663: PetscErrorCode EPSJDGetBOrth(EPS eps,PetscBool *borth)664: {
670: PetscTryMethod(eps,"EPSJDGetBOrth_C",(EPS,PetscBool*),(eps,borth));
671: return(0);
672: }
676: PETSC_EXTERN PetscErrorCode EPSCreate_JD(EPS eps)677: {
681: /* Load the Davidson solver */
682: EPSCreate_XD(eps);
683: EPSXDSetMethod(eps,DVD_METH_JD);
685: /* Overload the JD properties */
686: eps->ops->setfromoptions = EPSSetFromOptions_JD;
687: eps->ops->setup = EPSSetUp_JD;
688: eps->ops->destroy = EPSDestroy_JD;
690: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",EPSXDSetKrylovStart_XD);
691: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",EPSXDGetKrylovStart_XD);
692: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",EPSXDSetBlockSize_XD);
693: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",EPSXDGetBlockSize_XD);
694: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",EPSXDSetRestart_XD);
695: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",EPSXDGetRestart_XD);
696: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",EPSXDSetInitialSize_XD);
697: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",EPSXDGetInitialSize_XD);
698: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",EPSJDSetFix_JD);
699: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",EPSXDGetFix_XD);
700: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",EPSJDSetConstCorrectionTol_JD);
701: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",EPSJDGetConstCorrectionTol_JD);
702: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetWindowSizes_C",EPSXDSetWindowSizes_XD);
703: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetWindowSizes_C",EPSXDGetWindowSizes_XD);
704: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",EPSXDSetBOrth_XD);
705: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",EPSXDGetBOrth_XD);
706: return(0);
707: }