Actual source code: jd.c
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", submitted, 2013.
22: Last update: Jul 2012
24: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
25: SLEPc - Scalable Library for Eigenvalue Problem Computations
26: Copyright (c) 2002-2013, Universitat Politecnica de Valencia, Spain
28: This file is part of SLEPc.
30: SLEPc is free software: you can redistribute it and/or modify it under the
31: terms of version 3 of the GNU Lesser General Public License as published by
32: the Free Software Foundation.
34: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
35: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
36: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
37: more details.
39: You should have received a copy of the GNU Lesser General Public License
40: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
41: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
42: */
44: #include <slepc-private/epsimpl.h> /*I "slepceps.h" I*/
45: #include <../src/eps/impls/davidson/common/davidson.h>
49: PetscErrorCode EPSSetFromOptions_JD(EPS eps)
50: {
52: PetscBool flg,op;
53: PetscInt opi,opi0;
54: PetscReal opf;
55: KSP ksp;
56: EPSOrthType orth;
57: const char *orth_list[3] = {"I","B","B_opt"};
60: PetscOptionsHead("EPS Jacobi-Davidson (JD) Options");
62: EPSJDGetKrylovStart(eps,&op);
63: PetscOptionsBool("-eps_jd_krylov_start","Start the searching subspace with a krylov basis","EPSJDSetKrylovStart",op,&op,&flg);
64: if (flg) { EPSJDSetKrylovStart(eps,op); }
66: EPSJDGetBlockSize(eps,&opi);
67: PetscOptionsInt("-eps_jd_blocksize","Number vectors add to the searching subspace","EPSJDSetBlockSize",opi,&opi,&flg);
68: if (flg) { EPSJDSetBlockSize(eps,opi); }
70: EPSJDGetRestart(eps,&opi,&opi0);
71: PetscOptionsInt("-eps_jd_minv","Set the size of the searching subspace after restarting","EPSJDSetRestart",opi,&opi,&flg);
72: if (flg) { EPSJDSetRestart(eps,opi,opi0); }
74: PetscOptionsInt("-eps_jd_plusk","Set the number of saved eigenvectors from the previous iteration when restarting","EPSJDSetRestart",opi0,&opi0,&flg);
75: if (flg) { EPSJDSetRestart(eps,opi,opi0); }
77: EPSJDGetInitialSize(eps,&opi);
78: PetscOptionsInt("-eps_jd_initial_size","Set the initial size of the searching subspace","EPSJDSetInitialSize",opi,&opi,&flg);
79: if (flg) { EPSJDSetInitialSize(eps,opi); }
81: EPSJDGetFix(eps,&opf);
82: PetscOptionsReal("-eps_jd_fix","Set the tolerance for changing the target in the correction equation","EPSJDSetFix",opf,&opf,&flg);
83: if (flg) { EPSJDSetFix(eps,opf); }
85: EPSJDGetBOrth(eps,&orth);
86: PetscOptionsEList("-eps_jd_borth","orthogonalization used in the search subspace","EPSJDSetBOrth",orth_list,3,orth_list[orth-1],&opi,&flg);
87: if (flg) { EPSJDSetBOrth(eps,(EPSOrthType)(opi+1)); }
89: EPSJDGetConstCorrectionTol(eps,&op);
90: PetscOptionsBool("-eps_jd_const_correction_tol","Disable the dynamic stopping criterion when solving the correction equation","EPSJDSetConstCorrectionTol",op,&op,&flg);
91: if (flg) { EPSJDSetConstCorrectionTol(eps,op); }
93: EPSJDGetWindowSizes(eps,&opi,&opi0);
94: PetscOptionsInt("-eps_jd_pwindow","(Experimental!) Set the number of converged vectors in the projector","EPSJDSetWindowSizes",opi,&opi,&flg);
95: if (flg) { EPSJDSetWindowSizes(eps,opi,opi0); }
97: PetscOptionsInt("-eps_jd_qwindow","(Experimental!) Set the number of converged vectors in the projected problem","EPSJDSetWindowSizes",opi0,&opi0,&flg);
98: if (flg) { EPSJDSetWindowSizes(eps,opi,opi0); }
100: /* Set STPrecond as the default ST */
101: if (!((PetscObject)eps->st)->type_name) {
102: STSetType(eps->st,STPRECOND);
103: }
104: STPrecondSetKSPHasMat(eps->st,PETSC_FALSE);
106: /* Set the default options of the KSP */
107: STGetKSP(eps->st,&ksp);
108: if (!((PetscObject)ksp)->type_name) {
109: KSPSetType(ksp,KSPBCGSL);
110: KSPSetTolerances(ksp,1e-4,PETSC_DEFAULT,PETSC_DEFAULT,90);
111: }
112: PetscOptionsTail();
113: return(0);
114: }
118: PetscErrorCode EPSSetUp_JD(EPS eps)
119: {
121: PetscBool t;
122: KSP ksp;
125: /* Setup common for all davidson solvers */
126: EPSSetUp_XD(eps);
128: /* Set the default options of the KSP */
129: STGetKSP(eps->st,&ksp);
130: if (!((PetscObject)ksp)->type_name) {
131: KSPSetType(ksp,KSPBCGSL);
132: KSPSetTolerances(ksp,1e-4,PETSC_DEFAULT,PETSC_DEFAULT,90);
133: }
135: /* Check some constraints */
136: PetscObjectTypeCompare((PetscObject)ksp,KSPPREONLY,&t);
137: if (t) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"EPSJD does not work with KSPPREONLY");
138: return(0);
139: }
143: PetscErrorCode EPSDestroy_JD(EPS eps)
144: {
145: PetscErrorCode ierr;
148: PetscFree(eps->data);
149: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",NULL);
150: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",NULL);
151: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",NULL);
152: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",NULL);
153: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",NULL);
154: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",NULL);
155: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",NULL);
156: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",NULL);
157: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",NULL);
158: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",NULL);
159: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",NULL);
160: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",NULL);
161: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetWindowSizes_C",NULL);
162: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetWindowSizes_C",NULL);
163: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",NULL);
164: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",NULL);
165: return(0);
166: }
170: /*@
171: EPSJDSetKrylovStart - Activates or deactivates starting the searching
172: subspace with a Krylov basis.
174: Logically Collective on EPS
176: Input Parameters:
177: + eps - the eigenproblem solver context
178: - krylovstart - boolean flag
180: Options Database Key:
181: . -eps_jd_krylov_start - Activates starting the searching subspace with a
182: Krylov basis
184: Level: advanced
186: .seealso: EPSJDGetKrylovStart()
187: @*/
188: PetscErrorCode EPSJDSetKrylovStart(EPS eps,PetscBool krylovstart)
189: {
195: PetscTryMethod(eps,"EPSJDSetKrylovStart_C",(EPS,PetscBool),(eps,krylovstart));
196: return(0);
197: }
201: /*@
202: EPSJDGetKrylovStart - Returns a flag indicating if the searching subspace is started with a
203: Krylov basis.
205: Not Collective
207: Input Parameter:
208: . eps - the eigenproblem solver context
210: Output Parameters:
211: . krylovstart - boolean flag indicating if the searching subspace is started
212: with a Krylov basis
214: Level: advanced
216: .seealso: EPSJDGetKrylovStart()
217: @*/
218: PetscErrorCode EPSJDGetKrylovStart(EPS eps,PetscBool *krylovstart)
219: {
225: PetscTryMethod(eps,"EPSJDGetKrylovStart_C",(EPS,PetscBool*),(eps,krylovstart));
226: return(0);
227: }
231: /*@
232: EPSJDSetBlockSize - Sets the number of vectors to be added to the searching space
233: in every iteration.
235: Logically Collective on EPS
237: Input Parameters:
238: + eps - the eigenproblem solver context
239: - blocksize - number of vectors added to the search space in every iteration
241: Options Database Key:
242: . -eps_jd_blocksize - number of vectors added to the searching space every iteration
244: Level: advanced
246: .seealso: EPSJDSetKrylovStart()
247: @*/
248: PetscErrorCode EPSJDSetBlockSize(EPS eps,PetscInt blocksize)
249: {
255: PetscTryMethod(eps,"EPSJDSetBlockSize_C",(EPS,PetscInt),(eps,blocksize));
256: return(0);
257: }
261: /*@
262: EPSJDGetBlockSize - Returns the number of vectors to be added to the searching space
263: in every iteration.
265: Not Collective
267: Input Parameter:
268: . eps - the eigenproblem solver context
270: Output Parameter:
271: . blocksize - number of vectors added to the search space in every iteration
273: Level: advanced
275: .seealso: EPSJDSetBlockSize()
276: @*/
277: PetscErrorCode EPSJDGetBlockSize(EPS eps,PetscInt *blocksize)
278: {
284: PetscTryMethod(eps,"EPSJDGetBlockSize_C",(EPS,PetscInt*),(eps,blocksize));
285: return(0);
286: }
290: /*@
291: EPSJDGetRestart - Gets the number of vectors of the searching space after
292: restarting and the number of vectors saved from the previous iteration.
294: Not Collective
296: Input Parameter:
297: . eps - the eigenproblem solver context
299: Output Parameter:
300: + minv - number of vectors of the searching subspace after restarting
301: - plusk - number of vectors saved from the previous iteration
303: Level: advanced
305: .seealso: EPSJDSetRestart()
306: @*/
307: PetscErrorCode EPSJDGetRestart(EPS eps,PetscInt *minv,PetscInt *plusk)
308: {
313: PetscTryMethod(eps,"EPSJDGetRestart_C",(EPS,PetscInt*,PetscInt*),(eps,minv,plusk));
314: return(0);
315: }
319: /*@
320: EPSJDSetRestart - Sets the number of vectors of the searching space after
321: restarting and the number of vectors saved from the previous iteration.
323: Logically Collective on EPS
325: Input Parameters:
326: + eps - the eigenproblem solver context
327: . minv - number of vectors of the searching subspace after restarting
328: - plusk - number of vectors saved from the previous iteration
330: Options Database Keys:
331: + -eps_jd_minv - number of vectors of the searching subspace after restarting
332: - -eps_jd_plusk - number of vectors saved from the previous iteration
334: Level: advanced
336: .seealso: EPSJDGetRestart()
337: @*/
338: PetscErrorCode EPSJDSetRestart(EPS eps,PetscInt minv,PetscInt plusk)
339: {
346: PetscTryMethod(eps,"EPSJDSetRestart_C",(EPS,PetscInt,PetscInt),(eps,minv,plusk));
347: return(0);
348: }
352: /*@
353: EPSJDGetInitialSize - Returns the initial size of the searching space.
355: Not Collective
357: Input Parameter:
358: . eps - the eigenproblem solver context
360: Output Parameter:
361: . initialsize - number of vectors of the initial searching subspace
363: Notes:
364: If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
365: EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
366: provided vectors are not enough, the solver completes the subspace with
367: random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
368: gets the first vector provided by the user or, if not available, a random vector,
369: and expands the Krylov basis up to initialsize vectors.
371: Level: advanced
373: .seealso: EPSJDSetInitialSize(), EPSJDGetKrylovStart()
374: @*/
375: PetscErrorCode EPSJDGetInitialSize(EPS eps,PetscInt *initialsize)
376: {
382: PetscTryMethod(eps,"EPSJDGetInitialSize_C",(EPS,PetscInt*),(eps,initialsize));
383: return(0);
384: }
388: /*@
389: EPSJDSetInitialSize - Sets the initial size of the searching space.
391: Logically Collective on EPS
393: Input Parameters:
394: + eps - the eigenproblem solver context
395: - initialsize - number of vectors of the initial searching subspace
397: Options Database Key:
398: . -eps_jd_initial_size - number of vectors of the initial searching subspace
400: Notes:
401: If EPSJDGetKrylovStart() is PETSC_FALSE and the user provides vectors with
402: EPSSetInitialSpace(), up to initialsize vectors will be used; and if the
403: provided vectors are not enough, the solver completes the subspace with
404: random vectors. In the case of EPSJDGetKrylovStart() being PETSC_TRUE, the solver
405: gets the first vector provided by the user or, if not available, a random vector,
406: and expands the Krylov basis up to initialsize vectors.
408: Level: advanced
410: .seealso: EPSJDGetInitialSize(), EPSJDGetKrylovStart()
411: @*/
412: PetscErrorCode EPSJDSetInitialSize(EPS eps,PetscInt initialsize)
413: {
419: PetscTryMethod(eps,"EPSJDSetInitialSize_C",(EPS,PetscInt),(eps,initialsize));
420: return(0);
421: }
425: /*@
426: EPSJDGetFix - Returns the threshold for changing the target in the correction
427: equation.
429: Not Collective
431: Input Parameter:
432: . eps - the eigenproblem solver context
434: Output Parameter:
435: . fix - threshold for changing the target
437: Note:
438: The target in the correction equation is fixed at the first iterations.
439: When the norm of the residual vector is lower than the fix value,
440: the target is set to the corresponding eigenvalue.
442: Level: advanced
444: .seealso: EPSJDSetFix()
445: @*/
446: PetscErrorCode EPSJDGetFix(EPS eps,PetscReal *fix)
447: {
453: PetscTryMethod(eps,"EPSJDGetFix_C",(EPS,PetscReal*),(eps,fix));
454: return(0);
455: }
459: /*@
460: EPSJDSetFix - Sets the threshold for changing the target in the correction
461: equation.
463: Logically Collective on EPS
465: Input Parameters:
466: + eps - the eigenproblem solver context
467: - fix - threshold for changing the target
469: Options Database Key:
470: . -eps_jd_fix - the fix value
472: Note:
473: The target in the correction equation is fixed at the first iterations.
474: When the norm of the residual vector is lower than the fix value,
475: the target is set to the corresponding eigenvalue.
477: Level: advanced
479: .seealso: EPSJDGetFix()
480: @*/
481: PetscErrorCode EPSJDSetFix(EPS eps,PetscReal fix)
482: {
488: PetscTryMethod(eps,"EPSJDSetFix_C",(EPS,PetscReal),(eps,fix));
489: return(0);
490: }
494: /*@
495: EPSJDSetConstCorrectionTol - If true, deactivates the dynamic stopping criterion
496: (also called Newton) that sets the KSP relative tolerance
497: to 0.5**i, where i is the number of EPS iterations from the last converged value.
499: Logically Collective on EPS
501: Input Parameters:
502: + eps - the eigenproblem solver context
503: - constant - if false, the KSP relative tolerance is set to 0.5**i.
505: Options Database Key:
506: . -eps_jd_const_correction_tol - Deactivates the dynamic stopping criterion.
508: Level: advanced
510: .seealso: EPSJDGetConstCorrectionTol()
511: @*/
512: PetscErrorCode EPSJDSetConstCorrectionTol(EPS eps,PetscBool constant)
513: {
519: PetscTryMethod(eps,"EPSJDSetConstCorrectionTol_C",(EPS,PetscBool),(eps,constant));
520: return(0);
521: }
525: /*@
526: EPSJDGetConstCorrectionTol - Returns a flag indicating if the dynamic stopping is being used for
527: solving the correction equation. If the flag is false the KSP relative tolerance is set
528: to 0.5**i, where i is the number of EPS iterations from the last converged value.
530: Not Collective
532: Input Parameter:
533: . eps - the eigenproblem solver context
535: Output Parameters:
536: . constant - boolean flag indicating if the dynamic stopping criterion is not being used.
538: Level: advanced
540: .seealso: EPSJDGetConstCorrectionTol()
541: @*/
542: PetscErrorCode EPSJDGetConstCorrectionTol(EPS eps,PetscBool *constant)
543: {
549: PetscTryMethod(eps,"EPSJDGetConstCorrectionTol",(EPS,PetscBool*),(eps,constant));
550: return(0);
551: }
555: /*@
556: EPSJDGetWindowSizes - Gets the number of converged vectors in the projected
557: problem (or Rayleigh quotient) and in the projector employed in the correction
558: equation.
560: Not Collective
562: Input Parameter:
563: . eps - the eigenproblem solver context
565: Output Parameter:
566: + pwindow - number of converged vectors in the projector
567: - qwindow - number of converged vectors in the projected problem
569: Level: advanced
571: .seealso: EPSJDSetWindowSizes()
572: @*/
573: PetscErrorCode EPSJDGetWindowSizes(EPS eps,PetscInt *pwindow,PetscInt *qwindow)
574: {
579: PetscTryMethod(eps,"EPSJDGetWindowSizes_C",(EPS,PetscInt*,PetscInt*),(eps,pwindow,qwindow));
580: return(0);
581: }
585: /*@
586: EPSJDSetWindowSizes - Sets the number of converged vectors in the projected
587: problem (or Rayleigh quotient) and in the projector employed in the correction
588: equation.
590: Logically Collective on EPS
592: Input Parameters:
593: + eps - the eigenproblem solver context
594: . pwindow - number of converged vectors in the projector
595: - qwindow - number of converged vectors in the projected problem
597: Options Database Keys:
598: + -eps_jd_pwindow - set the number of converged vectors in the projector
599: - -eps_jd_qwindow - set the number of converged vectors in the projected problem
601: Level: advanced
603: .seealso: EPSJDGetWindowSizes()
604: @*/
605: PetscErrorCode EPSJDSetWindowSizes(EPS eps,PetscInt pwindow,PetscInt qwindow)
606: {
613: PetscTryMethod(eps,"EPSJDSetWindowSizes_C",(EPS,PetscInt,PetscInt),(eps,pwindow,qwindow));
614: return(0);
615: }
619: /*@
620: EPSJDSetBOrth - Selects the orthogonalization that will be used in the search
621: subspace in case of generalized Hermitian problems.
623: Logically Collective on EPS
625: Input Parameters:
626: + eps - the eigenproblem solver context
627: - borth - the kind of orthogonalization
629: Possible values:
630: The parameter 'borth' can have one of these values
632: + EPS_ORTH_I - orthogonalization of the search subspace
633: . EPS_ORTH_B - B-orthogonalization of the search subspace
634: - EPS_ORTH_BOPT - B-orthogonalization of the search subspace with an alternative method
636: Options Database Key:
637: . -eps_jd_borth - Set the orthogonalization used in the search subspace
639: Notes:
640: If borth is EPS_ORTH_B, the solver uses a variant of Gram-Schmidt (selected in
641: IP associated to the EPS) with the inner product defined by the matrix problem B.
642: If borth is EPS_ORTH_BOPT, it uses another variant of Gram-Schmidt that only performs
643: one matrix-vector product although more than one reorthogonalization would be done.
645: Level: advanced
647: .seealso: EPSJDGetBOrth()
648: @*/
649: PetscErrorCode EPSJDSetBOrth(EPS eps,EPSOrthType borth)
650: {
656: PetscTryMethod(eps,"EPSJDSetBOrth_C",(EPS,EPSOrthType),(eps,borth));
657: return(0);
658: }
662: /*@
663: EPSJDGetBOrth - Returns the orthogonalization used in the search
664: subspace in case of generalized Hermitian problems.
666: Not Collective
668: Input Parameter:
669: . eps - the eigenproblem solver context
671: Output Parameters:
672: . borth - the kind of orthogonalization
674: Notes:
675: See EPSJDSetBOrth() for possible values of 'borth'.
677: Level: advanced
679: .seealso: EPSJDSetBOrth(), EPSOrthType
680: @*/
681: PetscErrorCode EPSJDGetBOrth(EPS eps,EPSOrthType *borth)
682: {
688: PetscTryMethod(eps,"EPSJDGetBOrth_C",(EPS,EPSOrthType*),(eps,borth));
689: return(0);
690: }
694: PETSC_EXTERN PetscErrorCode EPSCreate_JD(EPS eps)
695: {
699: /* Load the Davidson solver */
700: EPSCreate_XD(eps);
701: EPSXDSetMethod(eps,DVD_METH_JD);
703: /* Overload the JD properties */
704: eps->ops->setfromoptions = EPSSetFromOptions_JD;
705: eps->ops->setup = EPSSetUp_JD;
706: eps->ops->destroy = EPSDestroy_JD;
708: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetKrylovStart_C",EPSXDSetKrylovStart_XD);
709: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetKrylovStart_C",EPSXDGetKrylovStart_XD);
710: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBlockSize_C",EPSXDSetBlockSize_XD);
711: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBlockSize_C",EPSXDGetBlockSize_XD);
712: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetRestart_C",EPSXDSetRestart_XD);
713: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetRestart_C",EPSXDGetRestart_XD);
714: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetInitialSize_C",EPSXDSetInitialSize_XD);
715: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetInitialSize_C",EPSXDGetInitialSize_XD);
716: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetFix_C",EPSJDSetFix_JD);
717: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetFix_C",EPSXDGetFix_XD);
718: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetConstCorrectionTol_C",EPSJDSetConstCorrectionTol_JD);
719: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetConstCorrectionTol_C",EPSJDGetConstCorrectionTol_JD);
720: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetWindowSizes_C",EPSXDSetWindowSizes_XD);
721: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetWindowSizes_C",EPSXDGetWindowSizes_XD);
722: PetscObjectComposeFunction((PetscObject)eps,"EPSJDSetBOrth_C",EPSXDSetBOrth_XD);
723: PetscObjectComposeFunction((PetscObject)eps,"EPSJDGetBOrth_C",EPSXDGetBOrth_XD);
724: return(0);
725: }