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-2011, 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 <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,QEPDefaultConverged,PETSC_NULL);}
83: PetscOptionsBoolGroupEnd("-qep_convergence_absolute","Absolute error convergence test","QEPSetConvergenceTest",&flg);
84: if (flg) {QEPSetConvergenceTest(qep,QEPAbsoluteConverged,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: IPSetFromOptions(qep->ip);
163: PetscRandomSetFromOptions(qep->rand);
164: return(0);
165: }
169: /*@
170: QEPGetTolerances - Gets the tolerance and maximum iteration count used
171: by the QEP convergence tests.
173: Not Collective
175: Input Parameter:
176: . qep - the quadratic eigensolver context
177:
178: Output Parameters:
179: + tol - the convergence tolerance
180: - maxits - maximum number of iterations
182: Notes:
183: The user can specify PETSC_NULL for any parameter that is not needed.
185: Level: intermediate
187: .seealso: QEPSetTolerances()
188: @*/
189: PetscErrorCode QEPGetTolerances(QEP qep,PetscReal *tol,PetscInt *maxits)
190: {
193: if (tol) *tol = qep->tol;
194: if (maxits) *maxits = qep->max_it;
195: return(0);
196: }
200: /*@
201: QEPSetTolerances - Sets the tolerance and maximum iteration count used
202: by the QEP convergence tests.
204: Logically Collective on QEP
206: Input Parameters:
207: + qep - the quadratic eigensolver context
208: . tol - the convergence tolerance
209: - maxits - maximum number of iterations to use
211: Options Database Keys:
212: + -qep_tol <tol> - Sets the convergence tolerance
213: - -qep_max_it <maxits> - Sets the maximum number of iterations allowed
215: Notes:
216: Use PETSC_IGNORE for an argument that need not be changed.
218: Use PETSC_DECIDE for maxits to assign a reasonably good value, which is
219: dependent on the solution method.
221: Level: intermediate
223: .seealso: QEPGetTolerances()
224: @*/
225: PetscErrorCode QEPSetTolerances(QEP qep,PetscReal tol,PetscInt maxits)
226: {
231: if (tol != PETSC_IGNORE) {
232: if (tol == PETSC_DEFAULT) {
233: qep->tol = PETSC_DEFAULT;
234: } else {
235: if (tol < 0.0) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
236: qep->tol = tol;
237: }
238: }
239: if (maxits != PETSC_IGNORE) {
240: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
241: qep->max_it = 0;
242: qep->setupcalled = 0;
243: } else {
244: if (maxits < 0) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
245: qep->max_it = maxits;
246: }
247: }
248: return(0);
249: }
253: /*@
254: QEPGetDimensions - Gets the number of eigenvalues to compute
255: and the dimension of the subspace.
257: Not Collective
259: Input Parameter:
260: . qep - the quadratic eigensolver context
261:
262: Output Parameters:
263: + nev - number of eigenvalues to compute
264: . ncv - the maximum dimension of the subspace to be used by the solver
265: - mpd - the maximum dimension allowed for the projected problem
267: Notes:
268: The user can specify PETSC_NULL for any parameter that is not needed.
270: Level: intermediate
272: .seealso: QEPSetDimensions()
273: @*/
274: PetscErrorCode QEPGetDimensions(QEP qep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
275: {
278: if (nev) *nev = qep->nev;
279: if (ncv) *ncv = qep->ncv;
280: if (mpd) *mpd = qep->mpd;
281: return(0);
282: }
286: /*@
287: QEPSetDimensions - Sets the number of eigenvalues to compute
288: and the dimension of the subspace.
290: Logically Collective on QEP
292: Input Parameters:
293: + qep - the quadratic eigensolver context
294: . nev - number of eigenvalues to compute
295: . ncv - the maximum dimension of the subspace to be used by the solver
296: - mpd - the maximum dimension allowed for the projected problem
298: Options Database Keys:
299: + -qep_nev <nev> - Sets the number of eigenvalues
300: . -qep_ncv <ncv> - Sets the dimension of the subspace
301: - -qep_mpd <mpd> - Sets the maximum projected dimension
303: Notes:
304: Use PETSC_IGNORE to retain the previous value of any parameter.
306: Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
307: dependent on the solution method.
309: The parameters ncv and mpd are intimately related, so that the user is advised
310: to set one of them at most. Normal usage is the following:
311: (a) In cases where nev is small, the user sets ncv (a reasonable default is 2*nev).
312: (b) In cases where nev is large, the user sets mpd.
314: The value of ncv should always be between nev and (nev+mpd), typically
315: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
316: a smaller value should be used.
318: Level: intermediate
320: .seealso: QEPGetDimensions()
321: @*/
322: PetscErrorCode QEPSetDimensions(QEP qep,PetscInt nev,PetscInt ncv,PetscInt mpd)
323: {
329: if (nev != PETSC_IGNORE) {
330: if (nev<1) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
331: qep->nev = nev;
332: qep->setupcalled = 0;
333: }
334: if (ncv != PETSC_IGNORE) {
335: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
336: qep->ncv = 0;
337: } else {
338: if (ncv<1) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
339: qep->ncv = ncv;
340: }
341: qep->setupcalled = 0;
342: }
343: if (mpd != PETSC_IGNORE) {
344: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
345: qep->mpd = 0;
346: } else {
347: if (mpd<1) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
348: qep->mpd = mpd;
349: }
350: }
351: return(0);
352: }
356: /*@
357: QEPSetWhichEigenpairs - Specifies which portion of the spectrum is
358: to be sought.
360: Logically Collective on QEP
362: Input Parameters:
363: + qep - eigensolver context obtained from QEPCreate()
364: - which - the portion of the spectrum to be sought
366: Possible values:
367: The parameter 'which' can have one of these values
368:
369: + QEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
370: . QEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
371: . QEP_LARGEST_REAL - largest real parts
372: . QEP_SMALLEST_REAL - smallest real parts
373: . QEP_LARGEST_IMAGINARY - largest imaginary parts
374: - QEP_SMALLEST_IMAGINARY - smallest imaginary parts
376: Options Database Keys:
377: + -qep_largest_magnitude - Sets largest eigenvalues in magnitude
378: . -qep_smallest_magnitude - Sets smallest eigenvalues in magnitude
379: . -qep_largest_real - Sets largest real parts
380: . -qep_smallest_real - Sets smallest real parts
381: . -qep_largest_imaginary - Sets largest imaginary parts
382: - -qep_smallest_imaginary - Sets smallest imaginary parts
384: Notes:
385: Not all eigensolvers implemented in QEP account for all the possible values
386: stated above. If SLEPc is compiled for real numbers QEP_LARGEST_IMAGINARY
387: and QEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
388: for eigenvalue selection.
389:
390: Level: intermediate
392: .seealso: QEPGetWhichEigenpairs(), QEPWhich
393: @*/
394: PetscErrorCode QEPSetWhichEigenpairs(QEP qep,QEPWhich which)
395: {
399: if (which!=PETSC_IGNORE) {
400: if (which==PETSC_DECIDE || which==PETSC_DEFAULT) qep->which = (QEPWhich)0;
401: else switch (which) {
402: case QEP_LARGEST_MAGNITUDE:
403: case QEP_SMALLEST_MAGNITUDE:
404: case QEP_LARGEST_REAL:
405: case QEP_SMALLEST_REAL:
406: case QEP_LARGEST_IMAGINARY:
407: case QEP_SMALLEST_IMAGINARY:
408: if (qep->which != which) {
409: qep->setupcalled = 0;
410: qep->which = which;
411: }
412: break;
413: default:
414: SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
415: }
416: }
417: return(0);
418: }
422: /*@C
423: QEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
424: sought.
426: Not Collective
428: Input Parameter:
429: . qep - eigensolver context obtained from QEPCreate()
431: Output Parameter:
432: . which - the portion of the spectrum to be sought
434: Notes:
435: See QEPSetWhichEigenpairs() for possible values of 'which'.
437: Level: intermediate
439: .seealso: QEPSetWhichEigenpairs(), QEPWhich
440: @*/
441: PetscErrorCode QEPGetWhichEigenpairs(QEP qep,QEPWhich *which)
442: {
446: *which = qep->which;
447: return(0);
448: }
452: /*@
453: QEPSetLeftVectorsWanted - Specifies which eigenvectors are required.
455: Logically Collective on QEP
457: Input Parameters:
458: + qep - the quadratic eigensolver context
459: - leftvecs - whether left eigenvectors are required or not
461: Options Database Keys:
462: . -qep_left_vectors <boolean> - Sets/resets the boolean flag 'leftvecs'
464: Notes:
465: If the user sets leftvecs=PETSC_TRUE then the solver uses a variant of
466: the algorithm that computes both right and left eigenvectors. This is
467: usually much more costly. This option is not available in all solvers.
469: Level: intermediate
471: .seealso: QEPGetLeftVectorsWanted(), QEPGetEigenvectorLeft()
472: @*/
473: PetscErrorCode QEPSetLeftVectorsWanted(QEP qep,PetscBool leftvecs)
474: {
478: if (qep->leftvecs != leftvecs) {
479: qep->leftvecs = leftvecs;
480: qep->setupcalled = 0;
481: }
482: return(0);
483: }
487: /*@C
488: QEPGetLeftVectorsWanted - Returns the flag indicating whether left
489: eigenvectors are required or not.
491: Not Collective
493: Input Parameter:
494: . qep - the eigensolver context
496: Output Parameter:
497: . leftvecs - the returned flag
499: Level: intermediate
501: .seealso: QEPSetLeftVectorsWanted()
502: @*/
503: PetscErrorCode QEPGetLeftVectorsWanted(QEP qep,PetscBool *leftvecs)
504: {
508: *leftvecs = qep->leftvecs;
509: return(0);
510: }
514: /*@
515: QEPGetScaleFactor - Gets the factor used for scaling the quadratic eigenproblem.
517: Not Collective
519: Input Parameter:
520: . qep - the quadratic eigensolver context
521:
522: Output Parameters:
523: . alpha - the scaling factor
525: Notes:
526: If the user did not specify a scaling factor, then after QEPSolve() the
527: default value is returned.
529: Level: intermediate
531: .seealso: QEPSetScaleFactor(), QEPSolve()
532: @*/
533: PetscErrorCode QEPGetScaleFactor(QEP qep,PetscReal *alpha)
534: {
538: *alpha = qep->sfactor;
539: return(0);
540: }
544: /*@
545: QEPSetScaleFactor - Sets the scaling factor to be used for scaling the
546: quadratic problem before attempting to solve.
548: Logically Collective on QEP
550: Input Parameters:
551: + qep - the quadratic eigensolver context
552: - alpha - the scaling factor
554: Options Database Keys:
555: . -qep_scale <alpha> - Sets the scaling factor
557: Notes:
558: For the problem (l^2*M + l*C + K)*x = 0, the effect of scaling is to work
559: with matrices (alpha^2*M, alpha*C, K), then scale the computed eigenvalue.
561: The default is to scale with alpha = norm(K)/norm(M).
563: Level: intermediate
565: .seealso: QEPGetScaleFactor()
566: @*/
567: PetscErrorCode QEPSetScaleFactor(QEP qep,PetscReal alpha)
568: {
572: if (alpha != PETSC_IGNORE) {
573: if (alpha == PETSC_DEFAULT || alpha == PETSC_DECIDE) {
574: qep->sfactor = 0.0;
575: } else {
576: if (alpha < 0.0) SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of alpha. Must be > 0");
577: qep->sfactor = alpha;
578: }
579: }
580: return(0);
581: }
585: /*@
586: QEPSetProblemType - Specifies the type of the quadratic eigenvalue problem.
588: Logically Collective on QEP
590: Input Parameters:
591: + qep - the quadratic eigensolver context
592: - type - a known type of quadratic eigenvalue problem
594: Options Database Keys:
595: + -qep_general - general problem with no particular structure
596: . -qep_hermitian - problem whose coefficient matrices are Hermitian
597: - -qep_gyroscopic - problem with Hamiltonian structure
598:
599: Notes:
600: Allowed values for the problem type are: general (QEP_GENERAL), Hermitian
601: (QEP_HERMITIAN), and gyroscopic (QEP_GYROSCOPIC).
603: This function is used to instruct SLEPc to exploit certain structure in
604: the quadratic eigenproblem. By default, no particular structure is assumed.
606: If the problem matrices are Hermitian (symmetric in the real case) or
607: Hermitian/skew-Hermitian then the solver can exploit this fact to perform
608: less operations or provide better stability.
610: Level: intermediate
612: .seealso: QEPSetOperators(), QEPSetType(), QEPGetProblemType(), QEPProblemType
613: @*/
614: PetscErrorCode QEPSetProblemType(QEP qep,QEPProblemType type)
615: {
619: if (type!=QEP_GENERAL && type!=QEP_HERMITIAN && type!=QEP_GYROSCOPIC)
620: SETERRQ(((PetscObject)qep)->comm,PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
621: qep->problem_type = type;
622: return(0);
623: }
627: /*@C
628: QEPGetProblemType - Gets the problem type from the QEP object.
630: Not Collective
632: Input Parameter:
633: . qep - the quadratic eigensolver context
635: Output Parameter:
636: . type - name of QEP problem type
638: Level: intermediate
640: .seealso: QEPSetProblemType(), QEPProblemType
641: @*/
642: PetscErrorCode QEPGetProblemType(QEP qep,QEPProblemType *type)
643: {
647: *type = qep->problem_type;
648: return(0);
649: }
653: /*@C
654: QEPSetConvergenceTest - Sets a function to compute the error estimate used in
655: the convergence test.
657: Logically Collective on QEP
659: Input Parameters:
660: + qep - eigensolver context obtained from QEPCreate()
661: . func - a pointer to the convergence test function
662: - ctx - a context pointer (the last parameter to the convergence test function)
664: Calling Sequence of func:
665: $ func(QEP qep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal* errest,void *ctx)
667: + qep - eigensolver context obtained from QEPCreate()
668: . eigr - real part of the eigenvalue
669: . eigi - imaginary part of the eigenvalue
670: . res - residual norm associated to the eigenpair
671: . errest - (output) computed error estimate
672: - ctx - optional context, as set by QEPSetConvergenceTest()
674: Note:
675: If the error estimate returned by the convergence test function is less than
676: the tolerance, then the eigenvalue is accepted as converged.
678: Level: advanced
680: .seealso: QEPSetTolerances()
681: @*/
683: {
686: qep->conv_func = func;
687: qep->conv_ctx = ctx;
688: return(0);
689: }
693: /*@
694: QEPSetTrackAll - Specifies if the solver must compute the residual of all
695: approximate eigenpairs or not.
697: Logically Collective on QEP
699: Input Parameters:
700: + qep - the eigensolver context
701: - trackall - whether compute all residuals or not
703: Notes:
704: If the user sets trackall=PETSC_TRUE then the solver explicitly computes
705: the residual for each eigenpair approximation. Computing the residual is
706: usually an expensive operation and solvers commonly compute the associated
707: residual to the first unconverged eigenpair.
708: The options '-qep_monitor_all' and '-qep_monitor_draw_all' automatically
709: activates this option.
711: Level: intermediate
713: .seealso: QEPGetTrackAll()
714: @*/
715: PetscErrorCode QEPSetTrackAll(QEP qep,PetscBool trackall)
716: {
720: qep->trackall = trackall;
721: return(0);
722: }
726: /*@
727: QEPGetTrackAll - Returns the flag indicating whether all residual norms must
728: be computed or not.
730: Not Collective
732: Input Parameter:
733: . qep - the eigensolver context
735: Output Parameter:
736: . trackall - the returned flag
738: Level: intermediate
740: .seealso: QEPSetTrackAll()
741: @*/
742: PetscErrorCode QEPGetTrackAll(QEP qep,PetscBool *trackall)
743: {
747: *trackall = qep->trackall;
748: return(0);
749: }
753: /*@C
754: QEPSetOptionsPrefix - Sets the prefix used for searching for all
755: QEP options in the database.
757: Logically Collective on QEP
759: Input Parameters:
760: + qep - the quadratic eigensolver context
761: - prefix - the prefix string to prepend to all QEP option requests
763: Notes:
764: A hyphen (-) must NOT be given at the beginning of the prefix name.
765: The first character of all runtime options is AUTOMATICALLY the
766: hyphen.
768: For example, to distinguish between the runtime options for two
769: different QEP contexts, one could call
770: .vb
771: QEPSetOptionsPrefix(qep1,"qeig1_")
772: QEPSetOptionsPrefix(qep2,"qeig2_")
773: .ve
775: Level: advanced
777: .seealso: QEPAppendOptionsPrefix(), QEPGetOptionsPrefix()
778: @*/
779: PetscErrorCode QEPSetOptionsPrefix(QEP qep,const char *prefix)
780: {
785: PetscObjectSetOptionsPrefix((PetscObject)qep,prefix);
786: return(0);
787: }
788:
791: /*@C
792: QEPAppendOptionsPrefix - Appends to the prefix used for searching for all
793: QEP options in the database.
795: Logically Collective on QEP
797: Input Parameters:
798: + qep - the quadratic eigensolver context
799: - prefix - the prefix string to prepend to all QEP option requests
801: Notes:
802: A hyphen (-) must NOT be given at the beginning of the prefix name.
803: The first character of all runtime options is AUTOMATICALLY the hyphen.
805: Level: advanced
807: .seealso: QEPSetOptionsPrefix(), QEPGetOptionsPrefix()
808: @*/
809: PetscErrorCode QEPAppendOptionsPrefix(QEP qep,const char *prefix)
810: {
815: PetscObjectAppendOptionsPrefix((PetscObject)qep,prefix);
816: return(0);
817: }
821: /*@C
822: QEPGetOptionsPrefix - Gets the prefix used for searching for all
823: QEP options in the database.
825: Not Collective
827: Input Parameters:
828: . qep - the quadratic eigensolver context
830: Output Parameters:
831: . prefix - pointer to the prefix string used is returned
833: Notes: On the fortran side, the user should pass in a string 'prefix' of
834: sufficient length to hold the prefix.
836: Level: advanced
838: .seealso: QEPSetOptionsPrefix(), QEPAppendOptionsPrefix()
839: @*/
840: PetscErrorCode QEPGetOptionsPrefix(QEP qep,const char *prefix[])
841: {
847: PetscObjectGetOptionsPrefix((PetscObject)qep,prefix);
848: return(0);
849: }