Actual source code: qepopts.c
1: /*
2: QEP 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-2012, Universitat Politecnica de Valencia, Spain
9: This file is part of SLEPc.
10:
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/qepimpl.h> /*I "slepcqep.h" I*/
29: /*@
30: QEPSetFromOptions - Sets QEP options from the options database.
31: This routine must be called before QEPSetUp() if the user is to be
32: allowed to set the solver type.
34: Collective on QEP
36: Input Parameters:
37: . qep - the quadratic eigensolver context
39: Notes:
40: To see all options, run your program with the -help option.
42: Level: beginner
43: @*/
44: PetscErrorCode QEPSetFromOptions(QEP qep)
45: {
46: PetscErrorCode ierr;
47: char type[256],monfilename[PETSC_MAX_PATH_LEN];
48: PetscBool flg,val;
49: PetscReal r;
50: PetscInt i,j,k;
51: PetscViewer monviewer;
52: SlepcConvMonitor ctx;
56: if (!QEPRegisterAllCalled) { QEPRegisterAll(PETSC_NULL); }
57: if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
58: PetscOptionsBegin(((PetscObject)qep)->comm,((PetscObject)qep)->prefix,"Quadratic Eigenvalue Problem (QEP) Solver Options","QEP");
59: PetscOptionsList("-qep_type","Quadratic Eigenvalue Problem method","QEPSetType",QEPList,(char*)(((PetscObject)qep)->type_name?((PetscObject)qep)->type_name:QEPLINEAR),type,256,&flg);
60: if (flg) {
61: QEPSetType(qep,type);
62: } else if (!((PetscObject)qep)->type_name) {
63: QEPSetType(qep,QEPLINEAR);
64: }
66: PetscOptionsBoolGroupBegin("-qep_general","general quadratic eigenvalue problem","QEPSetProblemType",&flg);
67: if (flg) {QEPSetProblemType(qep,QEP_GENERAL);}
68: PetscOptionsBoolGroup("-qep_hermitian","hermitian quadratic eigenvalue problem","QEPSetProblemType",&flg);
69: if (flg) {QEPSetProblemType(qep,QEP_HERMITIAN);}
70: PetscOptionsBoolGroupEnd("-qep_gyroscopic","gyroscopic quadratic eigenvalue problem","QEPSetProblemType",&flg);
71: if (flg) {QEPSetProblemType(qep,QEP_GYROSCOPIC);}
73: r = PETSC_IGNORE;
74: PetscOptionsReal("-qep_scale","Scale factor","QEPSetScaleFactor",qep->sfactor,&r,PETSC_NULL);
75: QEPSetScaleFactor(qep,r);
77: r = i = PETSC_IGNORE;
78: PetscOptionsInt("-qep_max_it","Maximum number of iterations","QEPSetTolerances",qep->max_it,&i,PETSC_NULL);
79: PetscOptionsReal("-qep_tol","Tolerance","QEPSetTolerances",qep->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:qep->tol,&r,PETSC_NULL);
80: QEPSetTolerances(qep,r,i);
81: PetscOptionsBoolGroupBegin("-qep_convergence_default","Default (relative error) convergence test","QEPSetConvergenceTest",&flg);
82: if (flg) {QEPSetConvergenceTest(qep,QEPConvergedDefault,PETSC_NULL);}
83: PetscOptionsBoolGroupEnd("-qep_convergence_absolute","Absolute error convergence test","QEPSetConvergenceTest",&flg);
84: if (flg) {QEPSetConvergenceTest(qep,QEPConvergedAbsolute,PETSC_NULL);}
86: i = j = k = PETSC_IGNORE;
87: PetscOptionsInt("-qep_nev","Number of eigenvalues to compute","QEPSetDimensions",qep->nev,&i,PETSC_NULL);
88: PetscOptionsInt("-qep_ncv","Number of basis vectors","QEPSetDimensions",qep->ncv,&j,PETSC_NULL);
89: PetscOptionsInt("-qep_mpd","Maximum dimension of projected problem","QEPSetDimensions",qep->mpd,&k,PETSC_NULL);
90: QEPSetDimensions(qep,i,j,k);
91:
92: /* -----------------------------------------------------------------------*/
93: /*
94: Cancels all monitors hardwired into code before call to QEPSetFromOptions()
95: */
96: flg = PETSC_FALSE;
97: PetscOptionsBool("-qep_monitor_cancel","Remove any hardwired monitor routines","QEPMonitorCancel",flg,&flg,PETSC_NULL);
98: if (flg) {
99: QEPMonitorCancel(qep);
100: }
101: /*
102: Prints approximate eigenvalues and error estimates at each iteration
103: */
104: PetscOptionsString("-qep_monitor","Monitor first unconverged approximate eigenvalue and error estimate","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
105: if (flg) {
106: PetscViewerASCIIOpen(((PetscObject)qep)->comm,monfilename,&monviewer);
107: QEPMonitorSet(qep,QEPMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
108: }
109: PetscOptionsString("-qep_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
110: if (flg) {
111: PetscNew(struct _n_SlepcConvMonitor,&ctx);
112: PetscViewerASCIIOpen(((PetscObject)qep)->comm,monfilename,&ctx->viewer);
113: QEPMonitorSet(qep,QEPMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
114: }
115: PetscOptionsString("-qep_monitor_all","Monitor approximate eigenvalues and error estimates","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
116: if (flg) {
117: PetscViewerASCIIOpen(((PetscObject)qep)->comm,monfilename,&monviewer);
118: QEPMonitorSet(qep,QEPMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
119: QEPSetTrackAll(qep,PETSC_TRUE);
120: }
121: flg = PETSC_FALSE;
122: PetscOptionsBool("-qep_monitor_draw","Monitor first unconverged approximate error estimate graphically","QEPMonitorSet",flg,&flg,PETSC_NULL);
123: if (flg) {
124: QEPMonitorSet(qep,QEPMonitorLG,PETSC_NULL,PETSC_NULL);
125: }
126: flg = PETSC_FALSE;
127: PetscOptionsBool("-qep_monitor_draw_all","Monitor error estimates graphically","QEPMonitorSet",flg,&flg,PETSC_NULL);
128: if (flg) {
129: QEPMonitorSet(qep,QEPMonitorLGAll,PETSC_NULL,PETSC_NULL);
130: QEPSetTrackAll(qep,PETSC_TRUE);
131: }
132: /* -----------------------------------------------------------------------*/
134: PetscOptionsBoolGroupBegin("-qep_largest_magnitude","compute largest eigenvalues in magnitude","QEPSetWhichEigenpairs",&flg);
135: if (flg) {QEPSetWhichEigenpairs(qep,QEP_LARGEST_MAGNITUDE);}
136: PetscOptionsBoolGroup("-qep_smallest_magnitude","compute smallest eigenvalues in magnitude","QEPSetWhichEigenpairs",&flg);
137: if (flg) {QEPSetWhichEigenpairs(qep,QEP_SMALLEST_MAGNITUDE);}
138: PetscOptionsBoolGroup("-qep_largest_real","compute largest real parts","QEPSetWhichEigenpairs",&flg);
139: if (flg) {QEPSetWhichEigenpairs(qep,QEP_LARGEST_REAL);}
140: PetscOptionsBoolGroup("-qep_smallest_real","compute smallest real parts","QEPSetWhichEigenpairs",&flg);
141: if (flg) {QEPSetWhichEigenpairs(qep,QEP_SMALLEST_REAL);}
142: PetscOptionsBoolGroup("-qep_largest_imaginary","compute largest imaginary parts","QEPSetWhichEigenpairs",&flg);
143: if (flg) {QEPSetWhichEigenpairs(qep,QEP_LARGEST_IMAGINARY);}
144: PetscOptionsBoolGroupEnd("-qep_smallest_imaginary","compute smallest imaginary parts","QEPSetWhichEigenpairs",&flg);
145: if (flg) {QEPSetWhichEigenpairs(qep,QEP_SMALLEST_IMAGINARY);}
147: PetscOptionsBool("-qep_left_vectors","Compute left eigenvectors also","QEPSetLeftVectorsWanted",qep->leftvecs,&val,&flg);
148: if (flg) {
149: QEPSetLeftVectorsWanted(qep,val);
150: }
152: PetscOptionsName("-qep_view","Print detailed information on solver used","QEPView",0);
153: PetscOptionsName("-qep_view_binary","Save the matrices associated to the eigenproblem","QEPSetFromOptions",0);
154: PetscOptionsName("-qep_plot_eigs","Make a plot of the computed eigenvalues","QEPSolve",0);
155:
156: if (qep->ops->setfromoptions) {
157: (*qep->ops->setfromoptions)(qep);
158: }
159: PetscObjectProcessOptionsHandlers((PetscObject)qep);
160: PetscOptionsEnd();
162: if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
163: IPSetFromOptions(qep->ip);
164: if (!qep->ds) { QEPGetDS(qep,&qep->ds); }
165: DSSetFromOptions(qep->ds);
166: PetscRandomSetFromOptions(qep->rand);
167: return(0);
168: }
172: /*@
173: QEPGetTolerances - Gets the tolerance and maximum iteration count used
174: by the QEP convergence tests.
176: Not Collective
178: Input Parameter:
179: . qep - the quadratic eigensolver context
180:
181: Output Parameters:
182: + tol - the convergence tolerance
183: - maxits - maximum number of iterations
185: Notes:
186: The user can specify PETSC_NULL for any parameter that is not needed.
188: Level: intermediate
190: .seealso: QEPSetTolerances()
191: @*/
192: PetscErrorCode QEPGetTolerances(QEP qep,PetscReal *tol,PetscInt *maxits)
193: {
196: if (tol) *tol = qep->tol;
197: if (maxits) *maxits = qep->max_it;
198: return(0);
199: }
203: /*@
204: QEPSetTolerances - Sets the tolerance and maximum iteration count used
205: by the QEP convergence tests.
207: Logically Collective on QEP
209: Input Parameters:
210: + qep - the quadratic eigensolver context
211: . tol - the convergence tolerance
212: - maxits - maximum number of iterations to use
214: Options Database Keys:
215: + -qep_tol <tol> - Sets the convergence tolerance
216: - -qep_max_it <maxits> - Sets the maximum number of iterations allowed
218: Notes:
219: Use PETSC_IGNORE for an argument that need not be changed.
221: Use PETSC_DECIDE for maxits to assign a reasonably good value, which is
222: dependent on the solution method.
224: Level: intermediate
226: .seealso: QEPGetTolerances()
227: @*/
228: PetscErrorCode QEPSetTolerances(QEP qep,PetscReal tol,PetscInt maxits)
229: {
234: if (tol != PETSC_IGNORE) {
235: if (tol == PETSC_DEFAULT) {
236: qep->tol = PETSC_DEFAULT;
237: } else {
238: if (tol < 0.0) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
239: qep->tol = tol;
240: }
241: }
242: if (maxits != PETSC_IGNORE) {
243: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
244: qep->max_it = 0;
245: qep->setupcalled = 0;
246: } else {
247: if (maxits < 0) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
248: qep->max_it = maxits;
249: }
250: }
251: return(0);
252: }
256: /*@
257: QEPGetDimensions - Gets the number of eigenvalues to compute
258: and the dimension of the subspace.
260: Not Collective
262: Input Parameter:
263: . qep - the quadratic eigensolver context
264:
265: Output Parameters:
266: + nev - number of eigenvalues to compute
267: . ncv - the maximum dimension of the subspace to be used by the solver
268: - mpd - the maximum dimension allowed for the projected problem
270: Notes:
271: The user can specify PETSC_NULL for any parameter that is not needed.
273: Level: intermediate
275: .seealso: QEPSetDimensions()
276: @*/
277: PetscErrorCode QEPGetDimensions(QEP qep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
278: {
281: if (nev) *nev = qep->nev;
282: if (ncv) *ncv = qep->ncv;
283: if (mpd) *mpd = qep->mpd;
284: return(0);
285: }
289: /*@
290: QEPSetDimensions - Sets the number of eigenvalues to compute
291: and the dimension of the subspace.
293: Logically Collective on QEP
295: Input Parameters:
296: + qep - the quadratic eigensolver context
297: . nev - number of eigenvalues to compute
298: . ncv - the maximum dimension of the subspace to be used by the solver
299: - mpd - the maximum dimension allowed for the projected problem
301: Options Database Keys:
302: + -qep_nev <nev> - Sets the number of eigenvalues
303: . -qep_ncv <ncv> - Sets the dimension of the subspace
304: - -qep_mpd <mpd> - Sets the maximum projected dimension
306: Notes:
307: Use PETSC_IGNORE to retain the previous value of any parameter.
309: Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
310: dependent on the solution method.
312: The parameters ncv and mpd are intimately related, so that the user is advised
313: to set one of them at most. Normal usage is that
314: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
315: (b) in cases where nev is large, the user sets mpd.
317: The value of ncv should always be between nev and (nev+mpd), typically
318: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
319: a smaller value should be used.
321: Level: intermediate
323: .seealso: QEPGetDimensions()
324: @*/
325: PetscErrorCode QEPSetDimensions(QEP qep,PetscInt nev,PetscInt ncv,PetscInt mpd)
326: {
332: if (nev != PETSC_IGNORE) {
333: if (nev<1) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
334: qep->nev = nev;
335: qep->setupcalled = 0;
336: }
337: if (ncv != PETSC_IGNORE) {
338: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
339: qep->ncv = 0;
340: } else {
341: if (ncv<1) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
342: qep->ncv = ncv;
343: }
344: qep->setupcalled = 0;
345: }
346: if (mpd != PETSC_IGNORE) {
347: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
348: qep->mpd = 0;
349: } else {
350: if (mpd<1) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
351: qep->mpd = mpd;
352: }
353: }
354: return(0);
355: }
359: /*@
360: QEPSetWhichEigenpairs - Specifies which portion of the spectrum is
361: to be sought.
363: Logically Collective on QEP
365: Input Parameters:
366: + qep - eigensolver context obtained from QEPCreate()
367: - which - the portion of the spectrum to be sought
369: Possible values:
370: The parameter 'which' can have one of these values
371:
372: + QEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
373: . QEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
374: . QEP_LARGEST_REAL - largest real parts
375: . QEP_SMALLEST_REAL - smallest real parts
376: . QEP_LARGEST_IMAGINARY - largest imaginary parts
377: - QEP_SMALLEST_IMAGINARY - smallest imaginary parts
379: Options Database Keys:
380: + -qep_largest_magnitude - Sets largest eigenvalues in magnitude
381: . -qep_smallest_magnitude - Sets smallest eigenvalues in magnitude
382: . -qep_largest_real - Sets largest real parts
383: . -qep_smallest_real - Sets smallest real parts
384: . -qep_largest_imaginary - Sets largest imaginary parts
385: - -qep_smallest_imaginary - Sets smallest imaginary parts
387: Notes:
388: Not all eigensolvers implemented in QEP account for all the possible values
389: stated above. If SLEPc is compiled for real numbers QEP_LARGEST_IMAGINARY
390: and QEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
391: for eigenvalue selection.
392:
393: Level: intermediate
395: .seealso: QEPGetWhichEigenpairs(), QEPWhich
396: @*/
397: PetscErrorCode QEPSetWhichEigenpairs(QEP qep,QEPWhich which)
398: {
402: if (which!=PETSC_IGNORE) {
403: if (which==PETSC_DECIDE || which==PETSC_DEFAULT) qep->which = (QEPWhich)0;
404: else switch (which) {
405: case QEP_LARGEST_MAGNITUDE:
406: case QEP_SMALLEST_MAGNITUDE:
407: case QEP_LARGEST_REAL:
408: case QEP_SMALLEST_REAL:
409: case QEP_LARGEST_IMAGINARY:
410: case QEP_SMALLEST_IMAGINARY:
411: if (qep->which != which) {
412: qep->setupcalled = 0;
413: qep->which = which;
414: }
415: break;
416: default:
417: SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
418: }
419: }
420: return(0);
421: }
425: /*@C
426: QEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
427: sought.
429: Not Collective
431: Input Parameter:
432: . qep - eigensolver context obtained from QEPCreate()
434: Output Parameter:
435: . which - the portion of the spectrum to be sought
437: Notes:
438: See QEPSetWhichEigenpairs() for possible values of 'which'.
440: Level: intermediate
442: .seealso: QEPSetWhichEigenpairs(), QEPWhich
443: @*/
444: PetscErrorCode QEPGetWhichEigenpairs(QEP qep,QEPWhich *which)
445: {
449: *which = qep->which;
450: return(0);
451: }
455: /*@
456: QEPSetLeftVectorsWanted - Specifies which eigenvectors are required.
458: Logically Collective on QEP
460: Input Parameters:
461: + qep - the quadratic eigensolver context
462: - leftvecs - whether left eigenvectors are required or not
464: Options Database Keys:
465: . -qep_left_vectors <boolean> - Sets/resets the boolean flag 'leftvecs'
467: Notes:
468: If the user sets leftvecs=PETSC_TRUE then the solver uses a variant of
469: the algorithm that computes both right and left eigenvectors. This is
470: usually much more costly. This option is not available in all solvers.
472: Level: intermediate
474: .seealso: QEPGetLeftVectorsWanted()
475: @*/
476: PetscErrorCode QEPSetLeftVectorsWanted(QEP qep,PetscBool leftvecs)
477: {
481: if (qep->leftvecs != leftvecs) {
482: qep->leftvecs = leftvecs;
483: qep->setupcalled = 0;
484: }
485: return(0);
486: }
490: /*@C
491: QEPGetLeftVectorsWanted - Returns the flag indicating whether left
492: eigenvectors are required or not.
494: Not Collective
496: Input Parameter:
497: . qep - the eigensolver context
499: Output Parameter:
500: . leftvecs - the returned flag
502: Level: intermediate
504: .seealso: QEPSetLeftVectorsWanted()
505: @*/
506: PetscErrorCode QEPGetLeftVectorsWanted(QEP qep,PetscBool *leftvecs)
507: {
511: *leftvecs = qep->leftvecs;
512: return(0);
513: }
517: /*@
518: QEPGetScaleFactor - Gets the factor used for scaling the quadratic eigenproblem.
520: Not Collective
522: Input Parameter:
523: . qep - the quadratic eigensolver context
524:
525: Output Parameters:
526: . alpha - the scaling factor
528: Notes:
529: If the user did not specify a scaling factor, then after QEPSolve() the
530: default value is returned.
532: Level: intermediate
534: .seealso: QEPSetScaleFactor(), QEPSolve()
535: @*/
536: PetscErrorCode QEPGetScaleFactor(QEP qep,PetscReal *alpha)
537: {
541: *alpha = qep->sfactor;
542: return(0);
543: }
547: /*@
548: QEPSetScaleFactor - Sets the scaling factor to be used for scaling the
549: quadratic problem before attempting to solve.
551: Logically Collective on QEP
553: Input Parameters:
554: + qep - the quadratic eigensolver context
555: - alpha - the scaling factor
557: Options Database Keys:
558: . -qep_scale <alpha> - Sets the scaling factor
560: Notes:
561: For the problem (l^2*M + l*C + K)*x = 0, the effect of scaling is to work
562: with matrices (alpha^2*M, alpha*C, K), then scale the computed eigenvalue.
564: The default is to scale with alpha = norm(K)/norm(M).
566: Level: intermediate
568: .seealso: QEPGetScaleFactor()
569: @*/
570: PetscErrorCode QEPSetScaleFactor(QEP qep,PetscReal alpha)
571: {
575: if (alpha != PETSC_IGNORE) {
576: if (alpha == PETSC_DEFAULT || alpha == PETSC_DECIDE) {
577: qep->sfactor = 0.0;
578: } else {
579: if (alpha < 0.0) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of alpha. Must be > 0");
580: qep->sfactor = alpha;
581: }
582: }
583: return(0);
584: }
588: /*@
589: QEPSetProblemType - Specifies the type of the quadratic eigenvalue problem.
591: Logically Collective on QEP
593: Input Parameters:
594: + qep - the quadratic eigensolver context
595: - type - a known type of quadratic eigenvalue problem
597: Options Database Keys:
598: + -qep_general - general problem with no particular structure
599: . -qep_hermitian - problem whose coefficient matrices are Hermitian
600: - -qep_gyroscopic - problem with Hamiltonian structure
601:
602: Notes:
603: Allowed values for the problem type are: general (QEP_GENERAL), Hermitian
604: (QEP_HERMITIAN), and gyroscopic (QEP_GYROSCOPIC).
606: This function is used to instruct SLEPc to exploit certain structure in
607: the quadratic eigenproblem. By default, no particular structure is assumed.
609: If the problem matrices are Hermitian (symmetric in the real case) or
610: Hermitian/skew-Hermitian then the solver can exploit this fact to perform
611: less operations or provide better stability.
613: Level: intermediate
615: .seealso: QEPSetOperators(), QEPSetType(), QEPGetProblemType(), QEPProblemType
616: @*/
617: PetscErrorCode QEPSetProblemType(QEP qep,QEPProblemType type)
618: {
622: if (type!=QEP_GENERAL && type!=QEP_HERMITIAN && type!=QEP_GYROSCOPIC)
623: SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
624: qep->problem_type = type;
625: return(0);
626: }
630: /*@C
631: QEPGetProblemType - Gets the problem type from the QEP object.
633: Not Collective
635: Input Parameter:
636: . qep - the quadratic eigensolver context
638: Output Parameter:
639: . type - name of QEP problem type
641: Level: intermediate
643: .seealso: QEPSetProblemType(), QEPProblemType
644: @*/
645: PetscErrorCode QEPGetProblemType(QEP qep,QEPProblemType *type)
646: {
650: *type = qep->problem_type;
651: return(0);
652: }
656: /*@C
657: QEPSetConvergenceTest - Sets a function to compute the error estimate used in
658: the convergence test.
660: Logically Collective on QEP
662: Input Parameters:
663: + qep - eigensolver context obtained from QEPCreate()
664: . func - a pointer to the convergence test function
665: - ctx - a context pointer (the last parameter to the convergence test function)
667: Calling Sequence of func:
668: $ func(QEP qep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal* errest,void *ctx)
670: + qep - eigensolver context obtained from QEPCreate()
671: . eigr - real part of the eigenvalue
672: . eigi - imaginary part of the eigenvalue
673: . res - residual norm associated to the eigenpair
674: . errest - (output) computed error estimate
675: - ctx - optional context, as set by QEPSetConvergenceTest()
677: Note:
678: If the error estimate returned by the convergence test function is less than
679: the tolerance, then the eigenvalue is accepted as converged.
681: Level: advanced
683: .seealso: QEPSetTolerances()
684: @*/
685: extern PetscErrorCode QEPSetConvergenceTest(QEP qep,PetscErrorCode (*func)(QEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx)
686: {
689: qep->conv_func = func;
690: qep->conv_ctx = ctx;
691: return(0);
692: }
696: /*@
697: QEPSetTrackAll - Specifies if the solver must compute the residual of all
698: approximate eigenpairs or not.
700: Logically Collective on QEP
702: Input Parameters:
703: + qep - the eigensolver context
704: - trackall - whether compute all residuals or not
706: Notes:
707: If the user sets trackall=PETSC_TRUE then the solver explicitly computes
708: the residual for each eigenpair approximation. Computing the residual is
709: usually an expensive operation and solvers commonly compute the associated
710: residual to the first unconverged eigenpair.
711: The options '-qep_monitor_all' and '-qep_monitor_draw_all' automatically
712: activates this option.
714: Level: intermediate
716: .seealso: QEPGetTrackAll()
717: @*/
718: PetscErrorCode QEPSetTrackAll(QEP qep,PetscBool trackall)
719: {
723: qep->trackall = trackall;
724: return(0);
725: }
729: /*@
730: QEPGetTrackAll - Returns the flag indicating whether all residual norms must
731: be computed or not.
733: Not Collective
735: Input Parameter:
736: . qep - the eigensolver context
738: Output Parameter:
739: . trackall - the returned flag
741: Level: intermediate
743: .seealso: QEPSetTrackAll()
744: @*/
745: PetscErrorCode QEPGetTrackAll(QEP qep,PetscBool *trackall)
746: {
750: *trackall = qep->trackall;
751: return(0);
752: }
756: /*@C
757: QEPSetOptionsPrefix - Sets the prefix used for searching for all
758: QEP options in the database.
760: Logically Collective on QEP
762: Input Parameters:
763: + qep - the quadratic eigensolver context
764: - prefix - the prefix string to prepend to all QEP option requests
766: Notes:
767: A hyphen (-) must NOT be given at the beginning of the prefix name.
768: The first character of all runtime options is AUTOMATICALLY the
769: hyphen.
771: For example, to distinguish between the runtime options for two
772: different QEP contexts, one could call
773: .vb
774: QEPSetOptionsPrefix(qep1,"qeig1_")
775: QEPSetOptionsPrefix(qep2,"qeig2_")
776: .ve
778: Level: advanced
780: .seealso: QEPAppendOptionsPrefix(), QEPGetOptionsPrefix()
781: @*/
782: PetscErrorCode QEPSetOptionsPrefix(QEP qep,const char *prefix)
783: {
788: if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
789: IPSetOptionsPrefix(qep->ip,prefix);
790: if (!qep->ds) { QEPGetDS(qep,&qep->ds); }
791: DSSetOptionsPrefix(qep->ds,prefix);
792: PetscObjectSetOptionsPrefix((PetscObject)qep,prefix);
793: return(0);
794: }
795:
798: /*@C
799: QEPAppendOptionsPrefix - Appends to the prefix used for searching for all
800: QEP options in the database.
802: Logically Collective on QEP
804: Input Parameters:
805: + qep - the quadratic eigensolver context
806: - prefix - the prefix string to prepend to all QEP option requests
808: Notes:
809: A hyphen (-) must NOT be given at the beginning of the prefix name.
810: The first character of all runtime options is AUTOMATICALLY the hyphen.
812: Level: advanced
814: .seealso: QEPSetOptionsPrefix(), QEPGetOptionsPrefix()
815: @*/
816: PetscErrorCode QEPAppendOptionsPrefix(QEP qep,const char *prefix)
817: {
819: PetscBool flg;
820: EPS eps;
824: if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
825: IPSetOptionsPrefix(qep->ip,prefix);
826: if (!qep->ds) { QEPGetDS(qep,&qep->ds); }
827: DSSetOptionsPrefix(qep->ds,prefix);
828: PetscObjectAppendOptionsPrefix((PetscObject)qep,prefix);
829: PetscObjectTypeCompare((PetscObject)qep,QEPLINEAR,&flg);
830: if (flg) {
831: QEPLinearGetEPS(qep,&eps);
832: EPSSetOptionsPrefix(eps,((PetscObject)qep)->prefix);
833: EPSAppendOptionsPrefix(eps,"qep_");
834: }
835: return(0);
836: }
840: /*@C
841: QEPGetOptionsPrefix - Gets the prefix used for searching for all
842: QEP options in the database.
844: Not Collective
846: Input Parameters:
847: . qep - the quadratic eigensolver context
849: Output Parameters:
850: . prefix - pointer to the prefix string used is returned
852: Notes: On the fortran side, the user should pass in a string 'prefix' of
853: sufficient length to hold the prefix.
855: Level: advanced
857: .seealso: QEPSetOptionsPrefix(), QEPAppendOptionsPrefix()
858: @*/
859: PetscErrorCode QEPGetOptionsPrefix(QEP qep,const char *prefix[])
860: {
866: PetscObjectGetOptionsPrefix((PetscObject)qep,prefix);
867: return(0);
868: }