Actual source code: jd.c
1: /*
2: SLEPc eigensolver: "jd"
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2011, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
9:
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include <private/epsimpl.h> /*I "slepceps.h" I*/
25: #include <private/stimpl.h> /*I "slepcst.h" I*/
26: #include <../src/eps/impls/davidson/common/davidson.h>
27: #include <slepcblaslapack.h>
29: PetscErrorCode EPSSetUp_JD(EPS eps);
30: PetscErrorCode EPSDestroy_JD(EPS eps);
35: PetscErrorCode EPSSetFromOptions_JD(EPS eps)
36: {
38: PetscBool flg,op;
39: PetscInt opi,opi0;
40: PetscReal opf;
41: KSP ksp;
44: PetscOptionsHead("EPS Jacobi-Davidson (JD) Options");
46: EPSJDGetKrylovStart(eps,&op);
47: PetscOptionsBool("-eps_jd_krylov_start","Start the searching subspace with a krylov basis","EPSJDSetKrylovStart",op,&op,&flg);
48: if(flg) { EPSJDSetKrylovStart(eps,op); }
49:
50: EPSJDGetBlockSize(eps,&opi);
51: PetscOptionsInt("-eps_jd_blocksize","Number vectors add to the searching subspace","EPSJDSetBlockSize",opi,&opi,&flg);
52: if(flg) { EPSJDSetBlockSize(eps,opi); }
54: EPSJDGetRestart(eps,&opi,&opi0);
55: PetscOptionsInt("-eps_jd_minv","Set the size of the searching subspace after restarting","EPSJDSetRestart",opi,&opi,&flg);
56: if(flg) { EPSJDSetRestart(eps,opi,opi0); }
58: PetscOptionsInt("-eps_jd_plusk","Set the number of saved eigenvectors from the previous iteration when restarting","EPSJDSetRestart",opi0,&opi0,&flg);
59: if(flg) { EPSJDSetRestart(eps,opi,opi0); }
61: EPSJDGetInitialSize(eps,&opi);
62: PetscOptionsInt("-eps_jd_initial_size","Set the initial size of the searching subspace","EPSJDSetInitialSize",opi,&opi,&flg);
63: if(flg) { EPSJDSetInitialSize(eps,opi); }
65: EPSJDGetFix(eps,&opf);
66: PetscOptionsReal("-eps_jd_fix","Set the tolerance for changing the target in the correction equation","EPSJDSetFix",opf,&opf,&flg);
67: if(flg) { EPSJDSetFix(eps,opf); }
69: EPSJDGetConstantCorrectionTolerance(eps,&op);
70: PetscOptionsBool("-eps_jd_constant_correction_tolerance","Disable the dynamic stopping criterion when solving the correction equation","EPSJDSetConstantCorrectionTolerance",op,&op,&flg);
71: if(flg) { EPSJDSetConstantCorrectionTolerance(eps,op); }
73: EPSJDGetWindowSizes(eps,&opi,&opi0);
74: PetscOptionsInt("-eps_jd_pwindow","(Experimental!) Set the number of converged vectors in the projector","EPSJDSetWindowSizes",opi,&opi,&flg);
75: if(flg) { EPSJDSetWindowSizes(eps,opi,opi0); }
77: PetscOptionsInt("-eps_jd_qwindow","(Experimental!) Set the number of converged vectors in the projected problem","EPSJDSetWindowSizes",opi0,&opi0,&flg);
78: if(flg) { EPSJDSetWindowSizes(eps,opi,opi0); }
80: PetscOptionsTail();
82: /* Set STPrecond as the default ST */
83: if (!((PetscObject)eps->OP)->type_name) {
84: STSetType(eps->OP,STPRECOND);
85: }
86: STPrecondSetKSPHasMat(eps->OP,PETSC_FALSE);
88: /* Set the default options of the KSP */
89: STGetKSP(eps->OP,&ksp);
90: if (!((PetscObject)ksp)->type_name) {
91: KSPSetType(ksp,KSPGMRES);
92: KSPSetTolerances(ksp,1e-4,PETSC_DEFAULT,PETSC_DEFAULT,90);
93: }
94:
95: return(0);
96: }
101: PetscErrorCode EPSSetUp_JD(EPS eps)
102: {
104: PetscBool t;
105: KSP ksp;
109: /* Setup common for all davidson solvers */
110: EPSSetUp_Davidson(eps);
112: /* Set the default options of the KSP */
113: STGetKSP(eps->OP,&ksp);
114: if (!((PetscObject)ksp)->type_name) {
115: KSPSetType(ksp,KSPGMRES);
116: KSPSetTolerances(ksp,1e-4,PETSC_DEFAULT,PETSC_DEFAULT,90);
117: }
118:
119: /* Check some constraints */
120: PetscTypeCompare((PetscObject)ksp,KSPPREONLY,&t);
121: if (t) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_SUP,"EPSJD does not work with KSPPREONLY");
122: return(0);
123: }
128: PetscErrorCode EPSCreate_JD(EPS eps)
129: {
133: /* Load the Davidson solver */
134: EPSCreate_Davidson(eps);
136: /* Overload the JD properties */
137: eps->ops->setfromoptions = EPSSetFromOptions_JD;
138: eps->ops->setup = EPSSetUp_JD;
139: eps->ops->destroy = EPSDestroy_JD;
141: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetKrylovStart_C","EPSDavidsonSetKrylovStart_Davidson",EPSDavidsonSetKrylovStart_Davidson);
142: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetKrylovStart_C","EPSDavidsonGetKrylovStart_Davidson",EPSDavidsonGetKrylovStart_Davidson);
143: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetBlockSize_C","EPSDavidsonSetBlockSize_Davidson",EPSDavidsonSetBlockSize_Davidson);
144: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetBlockSize_C","EPSDavidsonGetBlockSize_Davidson",EPSDavidsonGetBlockSize_Davidson);
145: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetRestart_C","EPSDavidsonSetRestart_Davidson",EPSDavidsonSetRestart_Davidson);
146: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetRestart_C","EPSDavidsonGetRestart_Davidson",EPSDavidsonGetRestart_Davidson);
147: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetInitialSize_C","EPSDavidsonSetInitialSize_Davidson",EPSDavidsonSetInitialSize_Davidson);
148: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetInitialSize_C","EPSDavidsonGetInitialSize_Davidson",EPSDavidsonGetInitialSize_Davidson);
149: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetFix_C","EPSDavidsonSetFix_Davidson",EPSDavidsonSetFix_Davidson);
150: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetFix_C","EPSDavidsonGetFix_Davidson",EPSDavidsonGetFix_Davidson);
151: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetConstantCorrectionTolerance_C","EPSDavidsonSetConstantCorrectionTolerance_Davidson",EPSDavidsonSetConstantCorrectionTolerance_Davidson);
152: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetConstantCorrectionTolerance_C","EPSDavidsonGetConstantCorrectionTolerance_Davidson",EPSDavidsonGetConstantCorrectionTolerance_Davidson);
153: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetWindowSizes_C","EPSDavidsonSetWindowSizes_Davidson",EPSDavidsonSetWindowSizes_Davidson);
154: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetWindowSizes_C","EPSDavidsonGetWindowSizes_Davidson",EPSDavidsonGetWindowSizes_Davidson);
155: return(0);
156: }
161: PetscErrorCode EPSDestroy_JD(EPS eps)
162: {
163: PetscErrorCode ierr;
166: PetscFree(eps->data);
167: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetKrylovStart_C","",PETSC_NULL);
168: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetKrylovStart_C","",PETSC_NULL);
169: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetBlockSize_C","",PETSC_NULL);
170: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetBlockSize_C","",PETSC_NULL);
171: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetRestart_C","",PETSC_NULL);
172: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetRestart_C","",PETSC_NULL);
173: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetInitialSize_C","",PETSC_NULL);
174: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetInitialSize_C","",PETSC_NULL);
175: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetFix_C","",PETSC_NULL);
176: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetFix_C","",PETSC_NULL);
177: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetConstantCorrectionTolerance_C","",PETSC_NULL);
178: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetConstantCorrectionTolerance_C","",PETSC_NULL);
179: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDSetWindowSizes_C","",PETSC_NULL);
180: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSJDGetWindowSizes_C","",PETSC_NULL);
181: return(0);
182: }
186: /*@
187: EPSJDSetKrylovStart - Activates or deactivates starting the searching
188: subspace with a Krylov basis.
190: Logically Collective on EPS
192: Input Parameters:
193: + eps - the eigenproblem solver context
194: - krylovstart - boolean flag
196: Options Database Key:
197: . -eps_jd_krylov_start - Activates starting the searching subspace with a
198: Krylov basis
199:
200: Level: advanced
202: .seealso: EPSJDGetKrylovStart()
203: @*/
204: PetscErrorCode EPSJDSetKrylovStart(EPS eps,PetscBool krylovstart)
205: {
211: PetscTryMethod(eps,"EPSJDSetKrylovStart_C",(EPS,PetscBool),(eps,krylovstart));
212: return(0);
213: }
217: /*@
218: EPSJDGetKrylovStart - Returns a flag indicating if the searching subspace is started with a
219: Krylov basis.
221: Not Collective
223: Input Parameter:
224: . eps - the eigenproblem solver context
226: Output Parameters:
227: . krylovstart - boolean flag indicating if the searching subspace is started
228: with a Krylov basis
230: Level: advanced
232: .seealso: EPSJDGetKrylovStart()
233: @*/
234: PetscErrorCode EPSJDGetKrylovStart(EPS eps,PetscBool *krylovstart)
235: {
241: PetscTryMethod(eps,"EPSJDGetKrylovStart_C",(EPS,PetscBool*),(eps,krylovstart));
242: return(0);
243: }
247: /*@
248: EPSJDSetBlockSize - Sets the number of vectors to be added to the searching space
249: in every iteration.
251: Logically Collective on EPS
253: Input Parameters:
254: + eps - the eigenproblem solver context
255: - blocksize - number of vectors added to the search space in every iteration
257: Options Database Key:
258: . -eps_jd_blocksize - number of vectors added to the searching space every iteration
259:
260: Level: advanced
262: .seealso: EPSJDSetKrylovStart()
263: @*/
264: PetscErrorCode EPSJDSetBlockSize(EPS eps,PetscInt blocksize)
265: {
271: PetscTryMethod(eps,"EPSJDSetBlockSize_C",(EPS,PetscInt),(eps,blocksize));
272: return(0);
273: }
277: /*@
278: EPSJDGetBlockSize - Returns the number of vectors to be added to the searching space
279: in every iteration.
281: Not Collective
283: Input Parameter:
284: . eps - the eigenproblem solver context
286: Output Parameter:
287: . blocksize - number of vectors added to the search space in every iteration
289: Level: advanced
291: .seealso: EPSJDSetBlockSize()
292: @*/
293: PetscErrorCode EPSJDGetBlockSize(EPS eps,PetscInt *blocksize)
294: {
300: PetscTryMethod(eps,"EPSJDGetBlockSize_C",(EPS,PetscInt*),(eps,blocksize));
301: return(0);
302: }
306: /*@
307: EPSJDGetRestart - Gets the number of vectors of the searching space after
308: restarting and the number of vectors saved from the previous iteration.
310: Not Collective
312: Input Parameter:
313: . eps - the eigenproblem solver context
315: Output Parameter:
316: + minv - number of vectors of the searching subspace after restarting
317: - plusk - number of vectors saved from the previous iteration
319: Level: advanced
321: .seealso: EPSJDSetRestart()
322: @*/
323: PetscErrorCode EPSJDGetRestart(EPS eps,PetscInt *minv,PetscInt *plusk)
324: {
331: PetscTryMethod(eps,"EPSJDGetRestart_C",(EPS,PetscInt*,PetscInt*),(eps,minv,plusk));
332: return(0);
333: }
337: /*@
338: EPSJDSetRestart - Sets the number of vectors of the searching space after
339: restarting and the number of vectors saved from the previous iteration.
341: Logically Collective on EPS
343: Input Parameters:
344: + eps - the eigenproblem solver context
345: . minv - number of vectors of the searching subspace after restarting
346: - plusk - number of vectors saved from the previous iteration
348: Options Database Keys:
349: + -eps_jd_minv - number of vectors of the searching subspace after restarting
350: - -eps_jd_plusk - number of vectors saved from the previous iteration
351:
352: Level: advanced
354: .seealso: EPSJDGetRestart()
355: @*/
356: PetscErrorCode EPSJDSetRestart(EPS eps,PetscInt minv,PetscInt plusk)
357: {
364: PetscTryMethod(eps,"EPSJDSetRestart_C",(EPS,PetscInt,PetscInt),(eps,minv,plusk));
365: return(0);
366: }
370: /*@
371: EPSJDGetInitialSize - Returns the initial size of the searching space.
373: Not Collective
375: Input Parameter:
376: . eps - the eigenproblem solver context
378: Output Parameter:
379: . initialsize - number of vectors of the initial searching subspace
381: Notes:
382: If EPSGDGetKrylovStart is PETSC_FALSE and the user provides vectors with
383: EPSSetInitialSpace, up to initialsize vectors will be used; and if the
384: provided vectors are not enough, the solver completes the subspace with
385: random vectors. In the case of EPSGDGetKrylovStart being PETSC_TRUE, the solver
386: gets the first vector provided by the user or, if not available, a random vector,
387: and expands the Krylov basis up to initialsize vectors.
389: Level: advanced
391: .seealso: EPSJDSetInitialSize(), EPSJDGetKrylovStart()
392: @*/
393: PetscErrorCode EPSJDGetInitialSize(EPS eps,PetscInt *initialsize)
394: {
400: PetscTryMethod(eps,"EPSJDGetInitialSize_C",(EPS,PetscInt*),(eps,initialsize));
401: return(0);
402: }
406: /*@
407: EPSJDSetInitialSize - Sets the initial size of the searching space.
409: Logically Collective on EPS
411: Input Parameters:
412: + eps - the eigenproblem solver context
413: - initialsize - number of vectors of the initial searching subspace
415: Options Database Key:
416: . -eps_jd_initial_size - number of vectors of the initial searching subspace
417:
418: Notes:
419: If EPSGDGetKrylovStart is PETSC_FALSE and the user provides vectors with
420: EPSSetInitialSpace, up to initialsize vectors will be used; and if the
421: provided vectors are not enough, the solver completes the subspace with
422: random vectors. In the case of EPSGDGetKrylovStart being PETSC_TRUE, the solver
423: gets the first vector provided by the user or, if not available, a random vector,
424: and expands the Krylov basis up to initialsize vectors.
426: Level: advanced
428: .seealso: EPSJDGetInitialSize(), EPSJDGetKrylovStart()
429: @*/
430: PetscErrorCode EPSJDSetInitialSize(EPS eps,PetscInt initialsize)
431: {
437: PetscTryMethod(eps,"EPSJDSetInitialSize_C",(EPS,PetscInt),(eps,initialsize));
438: return(0);
439: }
443: /*@
444: EPSJDGetFix - Returns the threshold for changing the target in the correction
445: equation.
447: Not Collective
449: Input Parameter:
450: . eps - the eigenproblem solver context
452: Output Parameter:
453: . fix - threshold for changing the target
455: Note:
456: The target in the correction equation is fixed at the first iterations.
457: When the norm of the residual vector is lower than the fix value,
458: the target is set to the corresponding eigenvalue.
460: Level: advanced
462: .seealso: EPSJDSetFix()
463: @*/
464: PetscErrorCode EPSJDGetFix(EPS eps,PetscReal *fix)
465: {
471: PetscTryMethod(eps,"EPSJDGetFix_C",(EPS,PetscReal*),(eps,fix));
472: return(0);
473: }
477: /*@
478: EPSJDSetFix - Sets the threshold for changing the target in the correction
479: equation.
481: Logically Collective on EPS
483: Input Parameters:
484: + eps - the eigenproblem solver context
485: - fix - threshold for changing the target
487: Options Database Key:
488: . -eps_jd_fix - the fix value
489:
490: Note:
491: The target in the correction equation is fixed at the first iterations.
492: When the norm of the residual vector is lower than the fix value,
493: the target is set to the corresponding eigenvalue.
495: Level: advanced
497: .seealso: EPSJDGetFix()
498: @*/
499: PetscErrorCode EPSJDSetFix(EPS eps,PetscReal fix)
500: {
506: PetscTryMethod(eps,"EPSJDSetFix_C",(EPS,PetscReal),(eps,fix));
507: return(0);
508: }
512: /*@
513: EPSJDSetConstantCorrectionTolerance - If true, deactivates the dynamic stopping criterion (also called Newton) that sets the KSP relative tolerance
514: to 0.5**i, where i is the number of EPS iterations from the last converged value.
516: Logically Collective on EPS
518: Input Parameters:
519: + eps - the eigenproblem solver context
520: - constant - if false, the KSP relative tolerance is set to 0.5**i.
522: Options Database Key:
523: . -eps_jd_constant_correction_tolerance - Deactivates the dynamic stopping criterion.
524:
525: Level: advanced
527: .seealso: EPSJDGetConstantCorrectionTolerance()
528: @*/
529: PetscErrorCode EPSJDSetConstantCorrectionTolerance(EPS eps,PetscBool constant)
530: {
536: PetscTryMethod(eps,"EPSJDSetConstantCorrectionTolerance_C",(EPS,PetscBool),(eps,constant));
537: return(0);
538: }
542: /*@
543: EPSJDGetConstantCorrectionTolerance - Returns a flag indicating if the dynamic stopping is being used for
544: solving the correction equation. If the flag is false the KSP relative tolerance is set
545: to 0.5**i, where i is the number of EPS iterations from the last converged value.
547: Not Collective
549: Input Parameter:
550: . eps - the eigenproblem solver context
552: Output Parameters:
553: . constant - boolean flag indicating if the dynamic stopping criterion is not being used.
555: Level: advanced
557: .seealso: EPSJDGetConstantCorrectionTolerance()
558: @*/
559: PetscErrorCode EPSJDGetConstantCorrectionTolerance(EPS eps,PetscBool *constant)
560: {
566: PetscTryMethod(eps,"EPSJDGetConstantCorrectionTolerance",(EPS,PetscBool*),(eps,constant));
567: return(0);
568: }
572: /*@
573: EPSJDGetWindowSizes - Gets the number of converged vectors in the projected
574: problem (or Rayleigh quotient) and in the projector employed in the correction
575: equation.
577: Not Collective
579: Input Parameter:
580: . eps - the eigenproblem solver context
582: Output Parameter:
583: + pwindow - number of converged vectors in the projector
584: - qwindow - number of converged vectors in the projected problem
586: Level: advanced
588: .seealso: EPSJDSetWindowSizes()
589: @*/
590: PetscErrorCode EPSJDGetWindowSizes(EPS eps,PetscInt *pwindow,PetscInt *qwindow)
591: {
598: PetscTryMethod(eps,"EPSJDGetWindowSizes_C",(EPS,PetscInt*,PetscInt*),(eps,pwindow,qwindow));
599: return(0);
600: }
604: /*@
605: EPSJDSetWindowSizes - Sets the number of converged vectors in the projected
606: problem (or Rayleigh quotient) and in the projector employed in the correction
607: equation.
609: Logically Collective on EPS
611: Input Parameters:
612: + eps - the eigenproblem solver context
613: . pwindow - number of converged vectors in the projector
614: - qwindow - number of converged vectors in the projected problem
616: Options Database Keys:
617: + -eps_jd_pwindow - set the number of converged vectors in the projector
618: - -eps_jd_qwindow - set the number of converged vectors in the projected problem
619:
620: Level: advanced
622: .seealso: EPSJDGetWindowSizes()
623: @*/
624: PetscErrorCode EPSJDSetWindowSizes(EPS eps,PetscInt pwindow,PetscInt qwindow)
625: {
632: PetscTryMethod(eps,"EPSJDSetWindowSizes_C",(EPS,PetscInt,PetscInt),(eps,pwindow,qwindow));
633: return(0);
634: }