1: /*
2: PEP routines related to options that can be set via the command-line
3: or procedurally.
5: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
6: SLEPc - Scalable Library for Eigenvalue Problem Computations
7: Copyright (c) 2002-2014, Universitat Politecnica de Valencia, Spain
9: This file is part of SLEPc.
11: SLEPc is free software: you can redistribute it and/or modify it under the
12: terms of version 3 of the GNU Lesser General Public License as published by
13: the Free Software Foundation.
15: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
16: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
18: more details.
20: You should have received a copy of the GNU Lesser General Public License
21: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
22: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
23: */
25: #include <slepc-private/pepimpl.h> /*I "slepcpep.h" I*/
29: /*@
30: PEPSetFromOptions - Sets PEP options from the options database.
31: This routine must be called before PEPSetUp() if the user is to be
32: allowed to set the solver type.
34: Collective on PEP 36: Input Parameters:
37: . pep - the polynomial eigensolver context
39: Notes:
40: To see all options, run your program with the -help option.
42: Level: beginner
43: @*/
44: PetscErrorCode PEPSetFromOptions(PEP pep) 45: {
46: PetscErrorCode ierr;
47: char type[256],monfilename[PETSC_MAX_PATH_LEN];
48: PetscBool flg,flg1,flg2,flg3,flg4;
49: PetscReal r,t;
50: PetscScalar s;
51: PetscInt i,j,k;
52: PetscViewer monviewer;
53: SlepcConvMonitor ctx;
57: if (!PEPRegisterAllCalled) { PEPRegisterAll(); }
58: PetscObjectOptionsBegin((PetscObject)pep);
59: PetscOptionsFList("-pep_type","Polynomial Eigenvalue Problem method","PEPSetType",PEPList,(char*)(((PetscObject)pep)->type_name?((PetscObject)pep)->type_name:PEPTOAR),type,256,&flg);
60: if (flg) {
61: PEPSetType(pep,type);
62: } else if (!((PetscObject)pep)->type_name) {
63: PEPSetType(pep,PEPTOAR);
64: }
66: PetscOptionsBoolGroupBegin("-pep_general","general polynomial eigenvalue problem","PEPSetProblemType",&flg);
67: if (flg) { PEPSetProblemType(pep,PEP_GENERAL); }
68: PetscOptionsBoolGroup("-pep_hermitian","hermitian polynomial eigenvalue problem","PEPSetProblemType",&flg);
69: if (flg) { PEPSetProblemType(pep,PEP_HERMITIAN); }
70: PetscOptionsBoolGroupEnd("-pep_gyroscopic","gyroscopic polynomial eigenvalue problem","PEPSetProblemType",&flg);
71: if (flg) { PEPSetProblemType(pep,PEP_GYROSCOPIC); }
73: PetscOptionsEnum("-pep_scale","Scaling strategy","PEPSetScale",PEPScaleTypes,(PetscEnum)pep->scale,(PetscEnum*)&pep->scale,NULL);
75: r = pep->sfactor;
76: PetscOptionsReal("-pep_scale_factor","Scale factor","PEPSetScale",pep->sfactor,&r,&flg1);
77: j = pep->sits;
78: PetscOptionsInt("-pep_scale_its","Number of iterations in diagonal scaling","PEPSetScale",pep->sits,&j,&flg2);
79: t = pep->slambda;
80: PetscOptionsReal("-pep_scale_lambda","Estimate of eigenvalue (modulus) for diagonal scaling","PEPSetScale",pep->slambda,&t,&flg3);
81: if (flg1 || flg2 || flg3) {
82: PEPSetScale(pep,pep->scale,r,j,t);
83: }
85: PetscOptionsEnum("-pep_extract","Extraction method","PEPSetExtract",PEPExtractTypes,(PetscEnum)pep->extract,(PetscEnum*)&pep->extract,NULL);
87: PetscOptionsEnum("-pep_refine","Iterative refinement method","PEPSetRefine",PEPRefineTypes,(PetscEnum)pep->refine,(PetscEnum*)&pep->refine,NULL);
89: i = pep->npart;
90: PetscOptionsInt("-pep_refine_partitions","Number of partitions of the communicator for iterative refinement","PEPSetRefine",pep->npart,&i,&flg1);
91: r = pep->rtol;
92: PetscOptionsReal("-pep_refine_tol","Tolerance for iterative refinement","PEPSetRefine",pep->rtol,&r,&flg2);
93: j = pep->rits;
94: PetscOptionsInt("-pep_refine_its","Maximum number of iterations for iterative refinement","PEPSetRefine",pep->rits,&j,&flg3);
95: flg = pep->schur;
96: PetscOptionsBool("-pep_refine_schur","Use Schur complement for iterative refinement","PEPSetRefine",pep->schur,&flg,&flg4);
97: if (flg1 || flg2 || flg3 || flg4) {
98: PEPSetRefine(pep,pep->refine,i,r,j,flg);
99: }
101: i = pep->max_it? pep->max_it: PETSC_DEFAULT;
102: PetscOptionsInt("-pep_max_it","Maximum number of iterations","PEPSetTolerances",pep->max_it,&i,&flg1);
103: r = pep->tol;
104: PetscOptionsReal("-pep_tol","Tolerance","PEPSetTolerances",pep->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:pep->tol,&r,&flg2);
105: if (flg1 || flg2) {
106: PEPSetTolerances(pep,r,i);
107: }
109: PetscOptionsBoolGroupBegin("-pep_conv_eig","Relative error convergence test","PEPSetConvergenceTest",&flg);
110: if (flg) { PEPSetConvergenceTest(pep,PEP_CONV_EIG); }
111: PetscOptionsBoolGroup("-pep_conv_norm","Convergence test relative to the eigenvalue and the matrix norms","PEPSetConvergenceTest",&flg);
112: if (flg) { PEPSetConvergenceTest(pep,PEP_CONV_NORM); }
113: PetscOptionsBoolGroup("-pep_conv_abs","Absolute error convergence test","PEPSetConvergenceTest",&flg);
114: if (flg) { PEPSetConvergenceTest(pep,PEP_CONV_ABS); }
115: PetscOptionsBoolGroupEnd("-pep_conv_user","User-defined convergence test","PEPSetConvergenceTest",&flg);
116: if (flg) { PEPSetConvergenceTest(pep,PEP_CONV_USER); }
118: i = pep->nev;
119: PetscOptionsInt("-pep_nev","Number of eigenvalues to compute","PEPSetDimensions",pep->nev,&i,&flg1);
120: j = pep->ncv? pep->ncv: PETSC_DEFAULT;
121: PetscOptionsInt("-pep_ncv","Number of basis vectors","PEPSetDimensions",pep->ncv,&j,&flg2);
122: k = pep->mpd? pep->mpd: PETSC_DEFAULT;
123: PetscOptionsInt("-pep_mpd","Maximum dimension of projected problem","PEPSetDimensions",pep->mpd,&k,&flg3);
124: if (flg1 || flg2 || flg3) {
125: PEPSetDimensions(pep,i,j,k);
126: }
128: PetscOptionsScalar("-pep_target","Value of the target","PEPSetTarget",pep->target,&s,&flg);
129: if (flg) {
130: PEPSetWhichEigenpairs(pep,PEP_TARGET_MAGNITUDE);
131: PEPSetTarget(pep,s);
132: }
134: PetscOptionsEnum("-pep_basis","Polynomial basis","PEPSetBasis",PEPBasisTypes,(PetscEnum)pep->basis,(PetscEnum*)&pep->basis,NULL);
136: /* -----------------------------------------------------------------------*/
137: /*
138: Cancels all monitors hardwired into code before call to PEPSetFromOptions()
139: */
140: flg = PETSC_FALSE;
141: PetscOptionsBool("-pep_monitor_cancel","Remove any hardwired monitor routines","PEPMonitorCancel",flg,&flg,NULL);
142: if (flg) {
143: PEPMonitorCancel(pep);
144: }
145: /*
146: Prints approximate eigenvalues and error estimates at each iteration
147: */
148: PetscOptionsString("-pep_monitor","Monitor first unconverged approximate eigenvalue and error estimate","PEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
149: if (flg) {
150: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)pep),monfilename,&monviewer);
151: PEPMonitorSet(pep,PEPMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
152: }
153: PetscOptionsString("-pep_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","PEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
154: if (flg) {
155: PetscNew(&ctx);
156: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)pep),monfilename,&ctx->viewer);
157: PEPMonitorSet(pep,PEPMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
158: }
159: PetscOptionsString("-pep_monitor_all","Monitor approximate eigenvalues and error estimates","PEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
160: if (flg) {
161: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)pep),monfilename,&monviewer);
162: PEPMonitorSet(pep,PEPMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
163: PEPSetTrackAll(pep,PETSC_TRUE);
164: }
165: flg = PETSC_FALSE;
166: PetscOptionsBool("-pep_monitor_lg","Monitor first unconverged approximate error estimate graphically","PEPMonitorSet",flg,&flg,NULL);
167: if (flg) {
168: PEPMonitorSet(pep,PEPMonitorLG,NULL,NULL);
169: }
170: flg = PETSC_FALSE;
171: PetscOptionsBool("-pep_monitor_lg_all","Monitor error estimates graphically","PEPMonitorSet",flg,&flg,NULL);
172: if (flg) {
173: PEPMonitorSet(pep,PEPMonitorLGAll,NULL,NULL);
174: PEPSetTrackAll(pep,PETSC_TRUE);
175: }
176: /* -----------------------------------------------------------------------*/
178: PetscOptionsBoolGroupBegin("-pep_largest_magnitude","compute largest eigenvalues in magnitude","PEPSetWhichEigenpairs",&flg);
179: if (flg) { PEPSetWhichEigenpairs(pep,PEP_LARGEST_MAGNITUDE); }
180: PetscOptionsBoolGroup("-pep_smallest_magnitude","compute smallest eigenvalues in magnitude","PEPSetWhichEigenpairs",&flg);
181: if (flg) { PEPSetWhichEigenpairs(pep,PEP_SMALLEST_MAGNITUDE); }
182: PetscOptionsBoolGroup("-pep_largest_real","compute largest real parts","PEPSetWhichEigenpairs",&flg);
183: if (flg) { PEPSetWhichEigenpairs(pep,PEP_LARGEST_REAL); }
184: PetscOptionsBoolGroup("-pep_smallest_real","compute smallest real parts","PEPSetWhichEigenpairs",&flg);
185: if (flg) { PEPSetWhichEigenpairs(pep,PEP_SMALLEST_REAL); }
186: PetscOptionsBoolGroup("-pep_largest_imaginary","compute largest imaginary parts","PEPSetWhichEigenpairs",&flg);
187: if (flg) { PEPSetWhichEigenpairs(pep,PEP_LARGEST_IMAGINARY); }
188: PetscOptionsBoolGroup("-pep_smallest_imaginary","compute smallest imaginary parts","PEPSetWhichEigenpairs",&flg);
189: if (flg) { PEPSetWhichEigenpairs(pep,PEP_SMALLEST_IMAGINARY); }
190: PetscOptionsBoolGroup("-pep_target_magnitude","compute nearest eigenvalues to target","PEPSetWhichEigenpairs",&flg);
191: if (flg) { PEPSetWhichEigenpairs(pep,PEP_TARGET_MAGNITUDE); }
192: PetscOptionsBoolGroup("-pep_target_real","compute eigenvalues with real parts close to target","PEPSetWhichEigenpairs",&flg);
193: if (flg) { PEPSetWhichEigenpairs(pep,PEP_TARGET_REAL); }
194: PetscOptionsBoolGroupEnd("-pep_target_imaginary","compute eigenvalues with imaginary parts close to target","PEPSetWhichEigenpairs",&flg);
195: if (flg) { PEPSetWhichEigenpairs(pep,PEP_TARGET_IMAGINARY); }
197: PetscOptionsName("-pep_view","Print detailed information on solver used","PEPView",0);
198: PetscOptionsName("-pep_plot_eigs","Make a plot of the computed eigenvalues","PEPSolve",0);
200: if (pep->ops->setfromoptions) {
201: (*pep->ops->setfromoptions)(pep);
202: }
203: PetscObjectProcessOptionsHandlers((PetscObject)pep);
204: PetscOptionsEnd();
206: if (!pep->V) { PEPGetBV(pep,&pep->V); }
207: BVSetFromOptions(pep->V);
208: if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
209: RGSetFromOptions(pep->rg);
210: if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
211: DSSetFromOptions(pep->ds);
212: if (!pep->st) { PEPGetST(pep,&pep->st); }
213: STSetFromOptions(pep->st);
214: PetscRandomSetFromOptions(pep->rand);
215: return(0);
216: }
220: /*@
221: PEPGetTolerances - Gets the tolerance and maximum iteration count used
222: by the PEP convergence tests.
224: Not Collective
226: Input Parameter:
227: . pep - the polynomial eigensolver context
229: Output Parameters:
230: + tol - the convergence tolerance
231: - maxits - maximum number of iterations
233: Notes:
234: The user can specify NULL for any parameter that is not needed.
236: Level: intermediate
238: .seealso: PEPSetTolerances()
239: @*/
240: PetscErrorCode PEPGetTolerances(PEP pep,PetscReal *tol,PetscInt *maxits)241: {
244: if (tol) *tol = pep->tol;
245: if (maxits) *maxits = pep->max_it;
246: return(0);
247: }
251: /*@
252: PEPSetTolerances - Sets the tolerance and maximum iteration count used
253: by the PEP convergence tests.
255: Logically Collective on PEP257: Input Parameters:
258: + pep - the polynomial eigensolver context
259: . tol - the convergence tolerance
260: - maxits - maximum number of iterations to use
262: Options Database Keys:
263: + -pep_tol <tol> - Sets the convergence tolerance
264: - -pep_max_it <maxits> - Sets the maximum number of iterations allowed
266: Notes:
267: Use PETSC_DEFAULT for either argument to assign a reasonably good value.
269: Level: intermediate
271: .seealso: PEPGetTolerances()
272: @*/
273: PetscErrorCode PEPSetTolerances(PEP pep,PetscReal tol,PetscInt maxits)274: {
279: if (tol == PETSC_DEFAULT) {
280: pep->tol = PETSC_DEFAULT;
281: pep->state = PEP_STATE_INITIAL;
282: } else {
283: if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
284: pep->tol = tol;
285: }
286: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
287: pep->max_it = 0;
288: pep->state = PEP_STATE_INITIAL;
289: } else {
290: if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
291: pep->max_it = maxits;
292: }
293: return(0);
294: }
298: /*@
299: PEPGetDimensions - Gets the number of eigenvalues to compute
300: and the dimension of the subspace.
302: Not Collective
304: Input Parameter:
305: . pep - the polynomial eigensolver context
307: Output Parameters:
308: + nev - number of eigenvalues to compute
309: . ncv - the maximum dimension of the subspace to be used by the solver
310: - mpd - the maximum dimension allowed for the projected problem
312: Notes:
313: The user can specify NULL for any parameter that is not needed.
315: Level: intermediate
317: .seealso: PEPSetDimensions()
318: @*/
319: PetscErrorCode PEPGetDimensions(PEP pep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)320: {
323: if (nev) *nev = pep->nev;
324: if (ncv) *ncv = pep->ncv;
325: if (mpd) *mpd = pep->mpd;
326: return(0);
327: }
331: /*@
332: PEPSetDimensions - Sets the number of eigenvalues to compute
333: and the dimension of the subspace.
335: Logically Collective on PEP337: Input Parameters:
338: + pep - the polynomial eigensolver context
339: . nev - number of eigenvalues to compute
340: . ncv - the maximum dimension of the subspace to be used by the solver
341: - mpd - the maximum dimension allowed for the projected problem
343: Options Database Keys:
344: + -pep_nev <nev> - Sets the number of eigenvalues
345: . -pep_ncv <ncv> - Sets the dimension of the subspace
346: - -pep_mpd <mpd> - Sets the maximum projected dimension
348: Notes:
349: Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
350: dependent on the solution method.
352: The parameters ncv and mpd are intimately related, so that the user is advised
353: to set one of them at most. Normal usage is that
354: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
355: (b) in cases where nev is large, the user sets mpd.
357: The value of ncv should always be between nev and (nev+mpd), typically
358: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
359: a smaller value should be used.
361: Level: intermediate
363: .seealso: PEPGetDimensions()
364: @*/
365: PetscErrorCode PEPSetDimensions(PEP pep,PetscInt nev,PetscInt ncv,PetscInt mpd)366: {
372: if (nev<1) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
373: pep->nev = nev;
374: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
375: pep->ncv = 0;
376: } else {
377: if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
378: pep->ncv = ncv;
379: }
380: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
381: pep->mpd = 0;
382: } else {
383: if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
384: pep->mpd = mpd;
385: }
386: pep->state = PEP_STATE_INITIAL;
387: return(0);
388: }
392: /*@
393: PEPSetWhichEigenpairs - Specifies which portion of the spectrum is
394: to be sought.
396: Logically Collective on PEP398: Input Parameters:
399: + pep - eigensolver context obtained from PEPCreate()
400: - which - the portion of the spectrum to be sought
402: Possible values:
403: The parameter 'which' can have one of these values
405: + PEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
406: . PEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
407: . PEP_LARGEST_REAL - largest real parts
408: . PEP_SMALLEST_REAL - smallest real parts
409: . PEP_LARGEST_IMAGINARY - largest imaginary parts
410: . PEP_SMALLEST_IMAGINARY - smallest imaginary parts
411: . PEP_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
412: . PEP_TARGET_REAL - eigenvalues with real part closest to target
413: - PEP_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
415: Options Database Keys:
416: + -pep_largest_magnitude - Sets largest eigenvalues in magnitude
417: . -pep_smallest_magnitude - Sets smallest eigenvalues in magnitude
418: . -pep_largest_real - Sets largest real parts
419: . -pep_smallest_real - Sets smallest real parts
420: . -pep_largest_imaginary - Sets largest imaginary parts
421: . -pep_smallest_imaginary - Sets smallest imaginary parts
422: . -pep_target_magnitude - Sets eigenvalues closest to target
423: . -pep_target_real - Sets real parts closest to target
424: - -pep_target_imaginary - Sets imaginary parts closest to target
426: Notes:
427: Not all eigensolvers implemented in PEP account for all the possible values
428: stated above. If SLEPc is compiled for real numbers PEP_LARGEST_IMAGINARY
429: and PEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
430: for eigenvalue selection.
432: Level: intermediate
434: .seealso: PEPGetWhichEigenpairs(), PEPWhich435: @*/
436: PetscErrorCode PEPSetWhichEigenpairs(PEP pep,PEPWhich which)437: {
441: switch (which) {
442: case PEP_LARGEST_MAGNITUDE:
443: case PEP_SMALLEST_MAGNITUDE:
444: case PEP_LARGEST_REAL:
445: case PEP_SMALLEST_REAL:
446: case PEP_LARGEST_IMAGINARY:
447: case PEP_SMALLEST_IMAGINARY:
448: case PEP_TARGET_MAGNITUDE:
449: case PEP_TARGET_REAL:
450: #if defined(PETSC_USE_COMPLEX)
451: case PEP_TARGET_IMAGINARY:
452: #endif
453: if (pep->which != which) {
454: pep->state = PEP_STATE_INITIAL;
455: pep->which = which;
456: }
457: break;
458: default:459: SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
460: }
461: return(0);
462: }
466: /*@C
467: PEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
468: sought.
470: Not Collective
472: Input Parameter:
473: . pep - eigensolver context obtained from PEPCreate()
475: Output Parameter:
476: . which - the portion of the spectrum to be sought
478: Notes:
479: See PEPSetWhichEigenpairs() for possible values of 'which'.
481: Level: intermediate
483: .seealso: PEPSetWhichEigenpairs(), PEPWhich484: @*/
485: PetscErrorCode PEPGetWhichEigenpairs(PEP pep,PEPWhich *which)486: {
490: *which = pep->which;
491: return(0);
492: }
496: /*@
497: PEPSetProblemType - Specifies the type of the polynomial eigenvalue problem.
499: Logically Collective on PEP501: Input Parameters:
502: + pep - the polynomial eigensolver context
503: - type - a known type of polynomial eigenvalue problem
505: Options Database Keys:
506: + -pep_general - general problem with no particular structure
507: . -pep_hermitian - problem whose coefficient matrices are Hermitian
508: - -pep_gyroscopic - problem with Hamiltonian structure
510: Notes:
511: Allowed values for the problem type are: general (PEP_GENERAL), Hermitian
512: (PEP_HERMITIAN), and gyroscopic (PEP_GYROSCOPIC).
514: This function is used to instruct SLEPc to exploit certain structure in
515: the polynomial eigenproblem. By default, no particular structure is assumed.
517: If the problem matrices are Hermitian (symmetric in the real case) or
518: Hermitian/skew-Hermitian then the solver can exploit this fact to perform
519: less operations or provide better stability.
521: Level: intermediate
523: .seealso: PEPSetOperators(), PEPSetType(), PEPGetProblemType(), PEPProblemType524: @*/
525: PetscErrorCode PEPSetProblemType(PEP pep,PEPProblemType type)526: {
530: if (type!=PEP_GENERAL && type!=PEP_HERMITIAN && type!=PEP_GYROSCOPIC) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
531: pep->problem_type = type;
532: return(0);
533: }
537: /*@C
538: PEPGetProblemType - Gets the problem type from the PEP object.
540: Not Collective
542: Input Parameter:
543: . pep - the polynomial eigensolver context
545: Output Parameter:
546: . type - name of PEP problem type
548: Level: intermediate
550: .seealso: PEPSetProblemType(), PEPProblemType551: @*/
552: PetscErrorCode PEPGetProblemType(PEP pep,PEPProblemType *type)553: {
557: *type = pep->problem_type;
558: return(0);
559: }
563: /*@
564: PEPSetBasis - Specifies the type of polynomial basis used to describe the
565: polynomial eigenvalue problem.
567: Logically Collective on PEP569: Input Parameters:
570: + pep - the polynomial eigensolver context
571: - basis - the type of polynomial basis
573: Options Database Key:
574: . -pep_basis <basis> - Select the basis type
576: Notes:
577: By default, the coefficient matrices passed via PEPSetOperators() are
578: expressed in the monomial basis, i.e.
579: P(lambda) = A_0 + lambda*A_1 + lambda^2*A_2 + ... + lambda^d*A_d.
580: Other polynomial bases may have better numerical behaviour, but the user
581: must then pass the coefficient matrices accordingly.
583: Level: intermediate
585: .seealso: PEPSetOperators(), PEPGetBasis(), PEPBasis586: @*/
587: PetscErrorCode PEPSetBasis(PEP pep,PEPBasis basis)588: {
592: pep->basis = basis;
593: return(0);
594: }
598: /*@C
599: PEPGetBasis - Gets the type of polynomial basis from the PEP object.
601: Not Collective
603: Input Parameter:
604: . pep - the polynomial eigensolver context
606: Output Parameter:
607: . basis - the polynomial basis
609: Level: intermediate
611: .seealso: PEPSetBasis(), PEPBasis612: @*/
613: PetscErrorCode PEPGetBasis(PEP pep,PEPBasis *basis)614: {
618: *basis = pep->basis;
619: return(0);
620: }
624: /*@
625: PEPSetTrackAll - Specifies if the solver must compute the residual of all
626: approximate eigenpairs or not.
628: Logically Collective on PEP630: Input Parameters:
631: + pep - the eigensolver context
632: - trackall - whether compute all residuals or not
634: Notes:
635: If the user sets trackall=PETSC_TRUE then the solver explicitly computes
636: the residual for each eigenpair approximation. Computing the residual is
637: usually an expensive operation and solvers commonly compute the associated
638: residual to the first unconverged eigenpair.
639: The options '-pep_monitor_all' and '-pep_monitor_lg_all' automatically
640: activate this option.
642: Level: intermediate
644: .seealso: PEPGetTrackAll()
645: @*/
646: PetscErrorCode PEPSetTrackAll(PEP pep,PetscBool trackall)647: {
651: pep->trackall = trackall;
652: return(0);
653: }
657: /*@
658: PEPGetTrackAll - Returns the flag indicating whether all residual norms must
659: be computed or not.
661: Not Collective
663: Input Parameter:
664: . pep - the eigensolver context
666: Output Parameter:
667: . trackall - the returned flag
669: Level: intermediate
671: .seealso: PEPSetTrackAll()
672: @*/
673: PetscErrorCode PEPGetTrackAll(PEP pep,PetscBool *trackall)674: {
678: *trackall = pep->trackall;
679: return(0);
680: }
684: /*@C
685: PEPSetConvergenceTestFunction - Sets a function to compute the error estimate
686: used in the convergence test.
688: Logically Collective on PEP690: Input Parameters:
691: + pep - eigensolver context obtained from PEPCreate()
692: . func - a pointer to the convergence test function
693: . ctx - [optional] context for private data for the convergence routine
694: - destroy - [optional] destructor for the context (may be NULL;
695: PETSC_NULL_FUNCTION in Fortran)
697: Calling Sequence of func:
698: $ func(PEP pep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
700: + pep - eigensolver context obtained from PEPCreate()
701: . eigr - real part of the eigenvalue
702: . eigi - imaginary part of the eigenvalue
703: . res - residual norm associated to the eigenpair
704: . errest - (output) computed error estimate
705: - ctx - optional context, as set by PEPSetConvergenceTest()
707: Note:
708: If the error estimate returned by the convergence test function is less than
709: the tolerance, then the eigenvalue is accepted as converged.
711: Level: advanced
713: .seealso: PEPSetConvergenceTest(), PEPSetTolerances()
714: @*/
715: PetscErrorCode PEPSetConvergenceTestFunction(PEP pep,PetscErrorCode (*func)(PEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))716: {
721: if (pep->convergeddestroy) {
722: (*pep->convergeddestroy)(pep->convergedctx);
723: }
724: pep->converged = func;
725: pep->convergeddestroy = destroy;
726: pep->convergedctx = ctx;
727: if (func == PEPConvergedEigRelative) pep->conv = PEP_CONV_EIG;
728: else if (func == PEPConvergedNormRelative) pep->conv = PEP_CONV_NORM;
729: else if (func == PEPConvergedAbsolute) pep->conv = PEP_CONV_ABS;
730: else pep->conv = PEP_CONV_USER;
731: return(0);
732: }
736: /*@
737: PEPSetConvergenceTest - Specifies how to compute the error estimate
738: used in the convergence test.
740: Logically Collective on PEP742: Input Parameters:
743: + pep - eigensolver context obtained from PEPCreate()
744: - conv - the type of convergence test
746: Options Database Keys:
747: + -pep_conv_abs - Sets the absolute convergence test
748: . -pep_conv_eig - Sets the convergence test relative to the eigenvalue
749: . -pep_conv_norm - Sets the convergence test relative to the matrix norms
750: - -pep_conv_user - Selects the user-defined convergence test
752: Note:
753: The parameter 'conv' can have one of these values
754: + PEP_CONV_ABS - absolute error ||r||
755: . PEP_CONV_EIG - error relative to the eigenvalue l, ||r||/|l|
756: . PEP_CONV_NORM - error relative to the matrix norms
757: - PEP_CONV_USER - function set by PEPSetConvergenceTestFunction()
759: Level: intermediate
761: .seealso: PEPGetConvergenceTest(), PEPSetConvergenceTestFunction(), PEPConv762: @*/
763: PetscErrorCode PEPSetConvergenceTest(PEP pep,PEPConv conv)764: {
768: switch (conv) {
769: case PEP_CONV_ABS: pep->converged = PEPConvergedAbsolute; break;
770: case PEP_CONV_EIG: pep->converged = PEPConvergedEigRelative; break;
771: case PEP_CONV_NORM: pep->converged = PEPConvergedNormRelative; break;
772: case PEP_CONV_USER: break;
773: default:774: SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
775: }
776: pep->conv = conv;
777: return(0);
778: }
782: /*@
783: PEPGetConvergenceTest - Gets the method used to compute the error estimate
784: used in the convergence test.
786: Not Collective
788: Input Parameters:
789: . pep - eigensolver context obtained from PEPCreate()
791: Output Parameters:
792: . conv - the type of convergence test
794: Level: intermediate
796: .seealso: PEPSetConvergenceTest(), PEPConv797: @*/
798: PetscErrorCode PEPGetConvergenceTest(PEP pep,PEPConv *conv)799: {
803: *conv = pep->conv;
804: return(0);
805: }
809: /*@
810: PEPSetScale - Specifies the scaling strategy to be used.
812: Logically Collective on PEP814: Input Parameters:
815: + pep - the eigensolver context
816: . scale - scaling strategy
817: . alpha - the scaling factor used in the scalar strategy
818: . its - number of iterations of the diagonal scaling algorithm
819: - lambda - approximation to wanted eigenvalues (modulus)
821: Options Database Keys:
822: + -pep_scale <type> - scaling type, one of <none,scalar,diagonal,both>
823: . -pep_scale_factor <alpha> - the scaling factor
824: . -pep_scale_its <its> - number of iterations
825: - -pep_scale_lambda <lambda> - approximation to eigenvalues
827: Notes:
828: There are two non-exclusive scaling strategies: scalar and diagonal.
830: In the scalar strategy, scaling is applied to the eigenvalue, that is,
831: mu = lambda/alpha is the new eigenvalue and all matrices are scaled
832: accordingly. After solving the scaled problem, the original lambda is
833: recovered. Parameter 'alpha' must be positive. Use PETSC_DECIDE to let
834: the solver compute a reasonable scaling factor.
836: In the diagonal strategy, the solver works implicitly with matrix Dr*A*Dl,
837: where Dr and Dl are appropriate diagonal matrices. This improves the accuracy
838: of the computed results in some cases. This option requires MATAIJ matrices.
839: The parameter 'its' is the number of iterations performed by the method.
840: Parameter 'lambda' must be positive. Use PETSC_DECIDE or set lambda = 1.0 if no
841: information about eigenvalues is available.
843: Level: intermediate
845: .seealso: PEPGetScale()
846: @*/
847: PetscErrorCode PEPSetScale(PEP pep,PEPScale scale,PetscReal alpha,PetscInt its,PetscReal lambda)848: {
852: pep->scale = scale;
853: if (scale==PEP_SCALE_SCALAR || scale==PEP_SCALE_BOTH) {
855: if (alpha == PETSC_DEFAULT || alpha == PETSC_DECIDE) {
856: pep->sfactor = 0.0;
857: pep->sfactor_set = PETSC_FALSE;
858: } else {
859: if (alpha<=0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of alpha. Must be > 0");
860: pep->sfactor = alpha;
861: pep->sfactor_set = PETSC_TRUE;
862: }
863: }
864: if (scale==PEP_SCALE_DIAGONAL || scale==PEP_SCALE_BOTH) {
867: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) pep->sits = 5;
868: else pep->sits = its;
869: if (lambda==PETSC_DECIDE || lambda==PETSC_DEFAULT) pep->slambda = 1.0;
870: else if (lambda<=0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of lambda. Must be > 0");
871: else pep->slambda = lambda;
872: }
873: pep->state = PEP_STATE_INITIAL;
874: return(0);
875: }
879: /*@
880: PEPGetScale - Gets the scaling strategy used by the PEP object, and the
881: associated parameters.
883: Not Collective
885: Input Parameter:
886: . pep - the eigensolver context
888: Output Parameters:
889: + scale - scaling strategy
890: . alpha - the scaling factor used in the scalar strategy
891: . its - number of iterations of the diagonal scaling algorithm
892: - lambda - approximation to wanted eigenvalues (modulus)
894: Level: intermediate
896: Note:
897: The user can specify NULL for any parameter that is not needed.
899: .seealso: PEPSetScale()
900: @*/
901: PetscErrorCode PEPGetScale(PEP pep,PEPScale *scale,PetscReal *alpha,PetscInt *its,PetscReal *lambda)902: {
905: if (scale) *scale = pep->scale;
906: if (alpha) *alpha = pep->sfactor;
907: if (its) *its = pep->sits;
908: if (lambda) *lambda = pep->slambda;
909: return(0);
910: }
914: /*@
915: PEPSetExtract - Specifies the extraction strategy to be used.
917: Logically Collective on PEP919: Input Parameters:
920: + pep - the eigensolver context
921: - extract - extraction strategy
923: Options Database Keys:
924: . -pep_extract <type> - extraction type, one of <norm,residual,structured>
926: Level: intermediate
928: .seealso: PEPGetExtract()
929: @*/
930: PetscErrorCode PEPSetExtract(PEP pep,PEPExtract extract)931: {
935: pep->extract = extract;
936: return(0);
937: }
941: /*@
942: PEPGetExtract - Gets the extraction strategy used by the PEP object.
944: Not Collective
946: Input Parameter:
947: . pep - the eigensolver context
949: Output Parameter:
950: . extract - extraction strategy
952: Level: intermediate
954: .seealso: PEPSetExtract()
955: @*/
956: PetscErrorCode PEPGetExtract(PEP pep,PEPExtract *extract)957: {
960: if (extract) *extract = pep->extract;
961: return(0);
962: }
966: /*@
967: PEPSetRefine - Specifies the refinement type (and options) to be used
968: after the solve.
970: Logically Collective on PEP972: Input Parameters:
973: + pep - the polynomial eigensolver context
974: . refine - refinement type
975: . npart - number of partitions of the communicator
976: . tol - the convergence tolerance
977: . its - maximum number of refinement iterations
978: - schur - boolean flag to activate the Schur complement approach
980: Options Database Keys:
981: + -pep_refine <type> - refinement type, one of <none,simple,multiple>
982: . -pep_refine_partitions <n> - the number of partitions
983: . -pep_refine_tol <tol> - the tolerance
984: . -pep_refine_its <its> - number of iterations
985: - -pep_refine_schur - to set the Schur complement approach
987: Notes:
988: By default, iterative refinement is disabled, since it may be very
989: costly. There are two possible refinement strategies: simple and multiple.
990: The simple approach performs iterative refinement on each of the
991: converged eigenpairs individually, whereas the multiple strategy works
992: with the invariant pair as a whole, refining all eigenpairs simultaneously.
993: The latter may be required for the case of multiple eigenvalues.
995: In some cases, especially when using direct solvers within the
996: iterative refinement method, it may be helpful for improved scalability
997: to split the communicator in several partitions. The npart parameter
998: indicates how many partitions to use (defaults to 1).
1000: The tol and its parameters specify the stopping criterion. In the simple
1001: method, refinement continues until the residual of each eigenpair is
1002: below the tolerance (tol defaults to the PEP tol, but may be set to a
1003: different value). In contrast, the multiple method simply performs its
1004: refinement iterations (just one by default).
1006: The schur flag is used to change the way in which linear systems are
1007: solved, so that a Schur complement approach is used instead of explicitly
1008: building the coefficient matrix.
1010: Level: intermediate
1012: .seealso: PEPGetRefine()
1013: @*/
1014: PetscErrorCode PEPSetRefine(PEP pep,PEPRefine refine,PetscInt npart,PetscReal tol,PetscInt its,PetscBool schur)1015: {
1017: PetscMPIInt size;
1026: pep->refine = refine;
1027: if (refine) { /* process parameters only if not REFINE_NONE */
1028: if (npart == PETSC_DEFAULT || npart == PETSC_DECIDE) {
1029: pep->npart = 1;
1030: } else {
1031: MPI_Comm_size(PetscObjectComm((PetscObject)pep),&size);
1032: if (npart<1 || npart>size) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of npart");
1033: pep->npart = npart;
1034: }
1035: if (tol == PETSC_DEFAULT || tol == PETSC_DECIDE) {
1036: pep->rtol = pep->tol;
1037: } else {
1038: if (tol<=0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
1039: pep->rtol = tol;
1040: }
1041: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) {
1042: pep->rits = PETSC_DEFAULT;
1043: } else {
1044: if (its<=0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be > 0");
1045: pep->rits = its;
1046: }
1047: pep->schur = schur;
1048: }
1049: pep->state = PEP_STATE_INITIAL;
1050: return(0);
1051: }
1055: /*@
1056: PEPGetRefine - Gets the refinement strategy used by the PEP object, and the
1057: associated parameters.
1059: Not Collective
1061: Input Parameter:
1062: . pep - the polynomial eigensolver context
1064: Output Parameters:
1065: + refine - refinement type
1066: . npart - number of partitions of the communicator
1067: . tol - the convergence tolerance
1068: . its - maximum number of refinement iterations
1069: - schur - whether the Schur complement approach is being used
1071: Level: intermediate
1073: Note:
1074: The user can specify NULL for any parameter that is not needed.
1076: .seealso: PEPSetRefine()
1077: @*/
1078: PetscErrorCode PEPGetRefine(PEP pep,PEPRefine *refine,PetscInt *npart,PetscReal *tol,PetscInt *its,PetscBool *schur)1079: {
1082: if (refine) *refine = pep->refine;
1083: if (npart) *npart = pep->npart;
1084: if (tol) *tol = pep->rtol;
1085: if (its) *its = pep->rits;
1086: if (schur) *schur = pep->schur;
1087: return(0);
1088: }
1092: /*@C
1093: PEPSetOptionsPrefix - Sets the prefix used for searching for all
1094: PEP options in the database.
1096: Logically Collective on PEP1098: Input Parameters:
1099: + pep - the polynomial eigensolver context
1100: - prefix - the prefix string to prepend to all PEP option requests
1102: Notes:
1103: A hyphen (-) must NOT be given at the beginning of the prefix name.
1104: The first character of all runtime options is AUTOMATICALLY the
1105: hyphen.
1107: For example, to distinguish between the runtime options for two
1108: different PEP contexts, one could call
1109: .vb
1110: PEPSetOptionsPrefix(pep1,"qeig1_")
1111: PEPSetOptionsPrefix(pep2,"qeig2_")
1112: .ve
1114: Level: advanced
1116: .seealso: PEPAppendOptionsPrefix(), PEPGetOptionsPrefix()
1117: @*/
1118: PetscErrorCode PEPSetOptionsPrefix(PEP pep,const char *prefix)1119: {
1124: if (!pep->st) { PEPGetST(pep,&pep->st); }
1125: STSetOptionsPrefix(pep->st,prefix);
1126: if (!pep->V) { PEPGetBV(pep,&pep->V); }
1127: BVSetOptionsPrefix(pep->V,prefix);
1128: if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
1129: DSSetOptionsPrefix(pep->ds,prefix);
1130: if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
1131: RGSetOptionsPrefix(pep->rg,prefix);
1132: PetscObjectSetOptionsPrefix((PetscObject)pep,prefix);
1133: return(0);
1134: }
1138: /*@C
1139: PEPAppendOptionsPrefix - Appends to the prefix used for searching for all
1140: PEP options in the database.
1142: Logically Collective on PEP1144: Input Parameters:
1145: + pep - the polynomial eigensolver context
1146: - prefix - the prefix string to prepend to all PEP option requests
1148: Notes:
1149: A hyphen (-) must NOT be given at the beginning of the prefix name.
1150: The first character of all runtime options is AUTOMATICALLY the hyphen.
1152: Level: advanced
1154: .seealso: PEPSetOptionsPrefix(), PEPGetOptionsPrefix()
1155: @*/
1156: PetscErrorCode PEPAppendOptionsPrefix(PEP pep,const char *prefix)1157: {
1159: PetscBool flg;
1160: EPS eps;
1164: if (!pep->st) { PEPGetST(pep,&pep->st); }
1165: STAppendOptionsPrefix(pep->st,prefix);
1166: if (!pep->V) { PEPGetBV(pep,&pep->V); }
1167: BVSetOptionsPrefix(pep->V,prefix);
1168: if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
1169: DSSetOptionsPrefix(pep->ds,prefix);
1170: if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
1171: RGSetOptionsPrefix(pep->rg,prefix);
1172: PetscObjectAppendOptionsPrefix((PetscObject)pep,prefix);
1173: PetscObjectTypeCompare((PetscObject)pep,PEPLINEAR,&flg);
1174: if (flg) {
1175: PEPLinearGetEPS(pep,&eps);
1176: EPSSetOptionsPrefix(eps,((PetscObject)pep)->prefix);
1177: EPSAppendOptionsPrefix(eps,"pep_");
1178: }
1179: return(0);
1180: }
1184: /*@C
1185: PEPGetOptionsPrefix - Gets the prefix used for searching for all
1186: PEP options in the database.
1188: Not Collective
1190: Input Parameters:
1191: . pep - the polynomial eigensolver context
1193: Output Parameters:
1194: . prefix - pointer to the prefix string used is returned
1196: Notes: On the fortran side, the user should pass in a string 'prefix' of
1197: sufficient length to hold the prefix.
1199: Level: advanced
1201: .seealso: PEPSetOptionsPrefix(), PEPAppendOptionsPrefix()
1202: @*/
1203: PetscErrorCode PEPGetOptionsPrefix(PEP pep,const char *prefix[])1204: {
1210: PetscObjectGetOptionsPrefix((PetscObject)pep,prefix);
1211: return(0);
1212: }