Actual source code: pjdopt.c
slepc-3.6.1 2015-09-03
1: /*
2: Options of polynomial JD solver.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
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 <slepc/private/pepimpl.h> /*I "slepcpep.h" I*/
25: #include pjdp.h
29: PetscErrorCode PEPJDSetTolerances_JD(PEP pep,PetscReal mtol,PetscReal htol,PetscReal stol)
30: {
31: PEP_JD *pjd = (PEP_JD*)pep->data;
34: if (mtol==PETSC_DEFAULT) pjd->mtol = 1e-5;
35: else {
36: if (mtol<0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mtol. Must be > 0");
37: pjd->mtol = mtol;
38: }
39: if (htol==PETSC_DEFAULT) pjd->htol = 1e-2;
40: else {
41: if (htol<0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of htol. Must be > 0");
42: pjd->htol = htol;
43: }
44: if (stol==PETSC_DEFAULT) pjd->stol = 1e-2;
45: else {
46: if (stol<0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of stol. Must be > 0");
47: pjd->stol = stol;
48: }
49: return(0);
50: }
54: /*@
55: PEPJDSetTolerances - Sets various tolerance parameters for the Jacobi-Davidson
56: method.
58: Logically Collective on PEP
60: Input Parameters:
61: + pep - the eigenproblem solver context
62: . mtol - the multiplicity tolerance
63: . htol - the tolerance for harmonic extraction
64: - stol - the tolerance for harmonic shift
66: Options Database Key:
67: + -pep_jd_mtol - Sets the multiplicity tolerance
68: . -pep_jd_htol - Sets the harmonic tolerance
69: - -pep_jd_stol - Sets the shift tolerance
71: Level: advanced
73: .seealso: PEPJDGetTolerances()
74: @*/
75: PetscErrorCode PEPJDSetTolerances(PEP pep,PetscReal mtol,PetscReal htol,PetscReal stol)
76: {
84: PetscTryMethod(pep,"PEPJDSetTolerances_C",(PEP,PetscReal,PetscReal,PetscReal),(pep,mtol,htol,stol));
85: return(0);
86: }
90: PetscErrorCode PEPJDGetTolerances_JD(PEP pep,PetscReal *mtol,PetscReal *htol,PetscReal *stol)
91: {
92: PEP_JD *pjd = (PEP_JD*)pep->data;
95: *mtol = pjd->mtol;
96: *htol = pjd->htol;
97: *stol = pjd->stol;
98: return(0);
99: }
103: /*@
104: PEPJDGetTolerances - Gets various tolerances in the Jacobi-Davidson method.
106: Not Collective
108: Input Parameter:
109: . pep - the eigenproblem solver context
111: Output Parameter:
112: + mtol - the multiplicity tolerance
113: . htol - the harmonic tolerance
114: - stol - the shift tolerance
116: Level: advanced
118: .seealso: PEPJDSetTolerances()
119: @*/
120: PetscErrorCode PEPJDGetTolerances(PEP pep,PetscReal *mtol,PetscReal *htol,PetscReal *stol)
121: {
129: PetscTryMethod(pep,"PEPJDGetTolerances_C",(PEP,PetscReal*,PetscReal*,PetscReal*),(pep,mtol,htol,stol));
130: return(0);
131: }
135: PetscErrorCode PEPJDSetRestart_JD(PEP pep,PetscReal keep)
136: {
137: PEP_JD *pjd = (PEP_JD*)pep->data;
140: if (keep==PETSC_DEFAULT) pjd->keep = 0.5;
141: else {
142: if (keep<0.1 || keep>0.9) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"The keep argument must be in the range [0.1,0.9]");
143: pjd->keep = keep;
144: }
145: return(0);
146: }
150: /*@
151: PEPJDSetRestart - Sets the restart parameter for the Jacobi-Davidson
152: method, in particular the proportion of basis vectors that must be kept
153: after restart.
155: Logically Collective on PEP
157: Input Parameters:
158: + pep - the eigenproblem solver context
159: - keep - the number of vectors to be kept at restart
161: Options Database Key:
162: . -pep_jd_restart - Sets the restart parameter
164: Notes:
165: Allowed values are in the range [0.1,0.9]. The default is 0.5.
167: Level: advanced
169: .seealso: PEPJDGetRestart()
170: @*/
171: PetscErrorCode PEPJDSetRestart(PEP pep,PetscReal keep)
172: {
178: PetscTryMethod(pep,"PEPJDSetRestart_C",(PEP,PetscReal),(pep,keep));
179: return(0);
180: }
184: PetscErrorCode PEPJDGetRestart_JD(PEP pep,PetscReal *keep)
185: {
186: PEP_JD *pjd = (PEP_JD*)pep->data;
189: *keep = pjd->keep;
190: return(0);
191: }
195: /*@
196: PEPJDGetRestart - Gets the restart parameter used in the Jacobi-Davidson method.
198: Not Collective
200: Input Parameter:
201: . pep - the eigenproblem solver context
203: Output Parameter:
204: . keep - the restart parameter
206: Level: advanced
208: .seealso: PEPJDSetRestart()
209: @*/
210: PetscErrorCode PEPJDGetRestart(PEP pep,PetscReal *keep)
211: {
217: PetscTryMethod(pep,"PEPJDGetRestart_C",(PEP,PetscReal*),(pep,keep));
218: return(0);
219: }
223: PetscErrorCode PEPSetFromOptions_JD(PetscOptions *PetscOptionsObject,PEP pep)
224: {
226: PEP_JD *pjd = (PEP_JD*)pep->data;
227: PetscBool flg,flg1,flg2,flg3;
228: PetscReal r1,r2,r3;
231: PetscOptionsHead(PetscOptionsObject,"PEP JD Options");
232: PetscOptionsReal("-pep_jd_restart","Proportion of vectors kept after restart","PEPJDSetRestart",0.5,&r1,&flg);
233: if (flg) {
234: PEPJDSetRestart(pep,r1);
235: }
237: r1 = pjd->mtol;
238: PetscOptionsReal("-pep_jd_mtol","Multiplicity tolerance","PEPJDSetTolerances",pjd->mtol,&r1,&flg1);
239: r2 = pjd->htol;
240: PetscOptionsReal("-pep_jd_htol","Harmonic tolerance","PEPJDSetTolerances",pjd->htol,&r2,&flg2);
241: r3 = pjd->stol;
242: PetscOptionsReal("-pep_jd_stol","Shift tolerance","PEPJDSetTolerances",pjd->stol,&r3,&flg3);
243: if (flg1 || flg2 || flg3) {
244: PEPJDSetTolerances(pep,r1,r2,r3);
245: }
246: PetscOptionsTail();
247: return(0);
248: }
252: PetscErrorCode PEPView_JD(PEP pep,PetscViewer viewer)
253: {
255: PEP_JD *pjd = (PEP_JD*)pep->data;
256: PetscBool isascii;
259: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
260: if (isascii) {
261: PetscViewerASCIIPrintf(viewer," JD: %d%% of basis vectors kept after restart\n",(int)(100*pjd->keep));
262: PetscViewerASCIIPrintf(viewer," JD: multiplicity tolerance = %g\n",(double)pjd->mtol);
263: PetscViewerASCIIPrintf(viewer," JD: harmonic tolerance = %g\n",(double)pjd->htol);
264: PetscViewerASCIIPrintf(viewer," JD: shift tolerance = %g\n",(double)pjd->stol);
265: }
266: return(0);
267: }