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-2010, Universidad 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
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: {
47: char type[256],monfilename[PETSC_MAX_PATH_LEN];
48: PetscTruth flg,val;
49: PetscReal r;
50: PetscInt i,j,k;
51: PetscViewerASCIIMonitor monviewer;
52: QEPMONITOR_CONV *ctx;
56: PetscOptionsBegin(((PetscObject)qep)->comm,((PetscObject)qep)->prefix,"Quadratic Eigenvalue Problem (QEP) Solver Options","QEP");
57: PetscOptionsList("-qep_type","Quadratic Eigenvalue Problem method","QEPSetType",QEPList,(char*)(((PetscObject)qep)->type_name?((PetscObject)qep)->type_name:QEPLINEAR),type,256,&flg);
58: if (flg) {
59: QEPSetType(qep,type);
60: } else if (!((PetscObject)qep)->type_name) {
61: QEPSetType(qep,QEPLINEAR);
62: }
64: PetscOptionsTruthGroupBegin("-qep_general","general quadratic eigenvalue problem","QEPSetProblemType",&flg);
65: if (flg) {QEPSetProblemType(qep,QEP_GENERAL);}
66: PetscOptionsTruthGroup("-qep_hermitian","hermitian quadratic eigenvalue problem","QEPSetProblemType",&flg);
67: if (flg) {QEPSetProblemType(qep,QEP_HERMITIAN);}
68: PetscOptionsTruthGroupEnd("-qep_gyroscopic","gyroscopic quadratic eigenvalue problem","QEPSetProblemType",&flg);
69: if (flg) {QEPSetProblemType(qep,QEP_GYROSCOPIC);}
71: r = PETSC_IGNORE;
72: PetscOptionsReal("-qep_scale","Scale factor","QEPSetScaleFactor",qep->sfactor,&r,PETSC_NULL);
73: QEPSetScaleFactor(qep,r);
75: r = i = PETSC_IGNORE;
76: PetscOptionsInt("-qep_max_it","Maximum number of iterations","QEPSetTolerances",qep->max_it,&i,PETSC_NULL);
77: PetscOptionsReal("-qep_tol","Tolerance","QEPSetTolerances",qep->tol,&r,PETSC_NULL);
78: QEPSetTolerances(qep,r,i);
79: PetscOptionsTruthGroupBegin("-qep_convergence_default","Default (relative error) convergence test","QEPSetConvergenceTest",&flg);
80: if (flg) {QEPSetConvergenceTest(qep,QEPDefaultConverged,PETSC_NULL);}
81: PetscOptionsTruthGroupEnd("-qep_convergence_absolute","Absolute error convergence test","QEPSetConvergenceTest",&flg);
82: if (flg) {QEPSetConvergenceTest(qep,QEPAbsoluteConverged,PETSC_NULL);}
84: i = j = k = PETSC_IGNORE;
85: PetscOptionsInt("-qep_nev","Number of eigenvalues to compute","QEPSetDimensions",qep->nev,&i,PETSC_NULL);
86: PetscOptionsInt("-qep_ncv","Number of basis vectors","QEPSetDimensions",qep->ncv,&j,PETSC_NULL);
87: PetscOptionsInt("-qep_mpd","Maximum dimension of projected problem","QEPSetDimensions",qep->mpd,&k,PETSC_NULL);
88: QEPSetDimensions(qep,i,j,k);
89:
90: /* -----------------------------------------------------------------------*/
91: /*
92: Cancels all monitors hardwired into code before call to QEPSetFromOptions()
93: */
94: flg = PETSC_FALSE;
95: PetscOptionsTruth("-qep_monitor_cancel","Remove any hardwired monitor routines","QEPMonitorCancel",flg,&flg,PETSC_NULL);
96: if (flg) {
97: QEPMonitorCancel(qep);
98: }
99: /*
100: Prints approximate eigenvalues and error estimates at each iteration
101: */
102: PetscOptionsString("-qep_monitor","Monitor first unconverged approximate eigenvalue and error estimate","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
103: if (flg) {
104: PetscViewerASCIIMonitorCreate(((PetscObject)qep)->comm,monfilename,((PetscObject)qep)->tablevel,&monviewer);
105: QEPMonitorSet(qep,QEPMonitorFirst,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
106: }
107: PetscOptionsString("-qep_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
108: if (flg) {
109: PetscNew(QEPMONITOR_CONV,&ctx);
110: PetscViewerASCIIMonitorCreate(((PetscObject)qep)->comm,monfilename,((PetscObject)qep)->tablevel,&ctx->viewer);
111: QEPMonitorSet(qep,QEPMonitorConverged,ctx,(PetscErrorCode (*)(void*))QEPMonitorDestroy_Converged);
112: }
113: PetscOptionsString("-qep_monitor_all","Monitor approximate eigenvalues and error estimates","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
114: if (flg) {
115: PetscViewerASCIIMonitorCreate(((PetscObject)qep)->comm,monfilename,((PetscObject)qep)->tablevel,&monviewer);
116: QEPMonitorSet(qep,QEPMonitorAll,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
117: QEPSetTrackAll(qep,PETSC_TRUE);
118: }
119: flg = PETSC_FALSE;
120: PetscOptionsTruth("-qep_monitor_draw","Monitor first unconverged approximate error estimate graphically","QEPMonitorSet",flg,&flg,PETSC_NULL);
121: if (flg) {
122: QEPMonitorSet(qep,QEPMonitorLG,PETSC_NULL,PETSC_NULL);
123: }
124: flg = PETSC_FALSE;
125: PetscOptionsTruth("-qep_monitor_draw_all","Monitor error estimates graphically","QEPMonitorSet",flg,&flg,PETSC_NULL);
126: if (flg) {
127: QEPMonitorSet(qep,QEPMonitorLGAll,PETSC_NULL,PETSC_NULL);
128: QEPSetTrackAll(qep,PETSC_TRUE);
129: }
130: /* -----------------------------------------------------------------------*/
132: PetscOptionsTruthGroupBegin("-qep_largest_magnitude","compute largest eigenvalues in magnitude","QEPSetWhichEigenpairs",&flg);
133: if (flg) {QEPSetWhichEigenpairs(qep,QEP_LARGEST_MAGNITUDE);}
134: PetscOptionsTruthGroup("-qep_smallest_magnitude","compute smallest eigenvalues in magnitude","QEPSetWhichEigenpairs",&flg);
135: if (flg) {QEPSetWhichEigenpairs(qep,QEP_SMALLEST_MAGNITUDE);}
136: PetscOptionsTruthGroup("-qep_largest_real","compute largest real parts","QEPSetWhichEigenpairs",&flg);
137: if (flg) {QEPSetWhichEigenpairs(qep,QEP_LARGEST_REAL);}
138: PetscOptionsTruthGroup("-qep_smallest_real","compute smallest real parts","QEPSetWhichEigenpairs",&flg);
139: if (flg) {QEPSetWhichEigenpairs(qep,QEP_SMALLEST_REAL);}
140: PetscOptionsTruthGroup("-qep_largest_imaginary","compute largest imaginary parts","QEPSetWhichEigenpairs",&flg);
141: if (flg) {QEPSetWhichEigenpairs(qep,QEP_LARGEST_IMAGINARY);}
142: PetscOptionsTruthGroupEnd("-qep_smallest_imaginary","compute smallest imaginary parts","QEPSetWhichEigenpairs",&flg);
143: if (flg) {QEPSetWhichEigenpairs(qep,QEP_SMALLEST_IMAGINARY);}
145: PetscOptionsTruth("-qep_left_vectors","Compute left eigenvectors also","QEPSetLeftVectorsWanted",qep->leftvecs,&val,&flg);
146: if (flg) {
147: QEPSetLeftVectorsWanted(qep,val);
148: }
150: PetscOptionsName("-qep_view","Print detailed information on solver used","QEPView",0);
151: PetscOptionsName("-qep_view_binary","Save the matrices associated to the eigenproblem","QEPSetFromOptions",0);
152: PetscOptionsName("-qep_plot_eigs","Make a plot of the computed eigenvalues","QEPSolve",0);
153:
154: if (qep->ops->setfromoptions) {
155: (*qep->ops->setfromoptions)(qep);
156: }
157: PetscOptionsEnd();
159: IPSetFromOptions(qep->ip);
160: return(0);
161: }
165: /*@
166: QEPGetTolerances - Gets the tolerance and maximum iteration count used
167: by the QEP convergence tests.
169: Not Collective
171: Input Parameter:
172: . qep - the quadratic eigensolver context
173:
174: Output Parameters:
175: + tol - the convergence tolerance
176: - maxits - maximum number of iterations
178: Notes:
179: The user can specify PETSC_NULL for any parameter that is not needed.
181: Level: intermediate
183: .seealso: QEPSetTolerances()
184: @*/
185: PetscErrorCode QEPGetTolerances(QEP qep,PetscReal *tol,PetscInt *maxits)
186: {
189: if (tol) *tol = qep->tol;
190: if (maxits) *maxits = qep->max_it;
191: return(0);
192: }
196: /*@
197: QEPSetTolerances - Sets the tolerance and maximum iteration count used
198: by the QEP convergence tests.
200: Collective on QEP
202: Input Parameters:
203: + qep - the quadratic eigensolver context
204: . tol - the convergence tolerance
205: - maxits - maximum number of iterations to use
207: Options Database Keys:
208: + -qep_tol <tol> - Sets the convergence tolerance
209: - -qep_max_it <maxits> - Sets the maximum number of iterations allowed
211: Notes:
212: Use PETSC_IGNORE for an argument that need not be changed.
214: Use PETSC_DECIDE for maxits to assign a reasonably good value, which is
215: dependent on the solution method.
217: Level: intermediate
219: .seealso: QEPGetTolerances()
220: @*/
221: PetscErrorCode QEPSetTolerances(QEP qep,PetscReal tol,PetscInt maxits)
222: {
225: if (tol != PETSC_IGNORE) {
226: if (tol == PETSC_DEFAULT) {
227: qep->tol = 1e-7;
228: } else {
229: if (tol < 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
230: qep->tol = tol;
231: }
232: }
233: if (maxits != PETSC_IGNORE) {
234: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
235: qep->max_it = 0;
236: qep->setupcalled = 0;
237: } else {
238: if (maxits < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
239: qep->max_it = maxits;
240: }
241: }
242: return(0);
243: }
247: /*@
248: QEPGetDimensions - Gets the number of eigenvalues to compute
249: and the dimension of the subspace.
251: Not Collective
253: Input Parameter:
254: . qep - the quadratic eigensolver context
255:
256: Output Parameters:
257: + nev - number of eigenvalues to compute
258: . ncv - the maximum dimension of the subspace to be used by the solver
259: - mpd - the maximum dimension allowed for the projected problem
261: Notes:
262: The user can specify PETSC_NULL for any parameter that is not needed.
264: Level: intermediate
266: .seealso: QEPSetDimensions()
267: @*/
268: PetscErrorCode QEPGetDimensions(QEP qep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
269: {
272: if (nev) *nev = qep->nev;
273: if (ncv) *ncv = qep->ncv;
274: if (mpd) *mpd = qep->mpd;
275: return(0);
276: }
280: /*@
281: QEPSetDimensions - Sets the number of eigenvalues to compute
282: and the dimension of the subspace.
284: Collective on QEP
286: Input Parameters:
287: + qep - the quadratic eigensolver context
288: . nev - number of eigenvalues to compute
289: . ncv - the maximum dimension of the subspace to be used by the solver
290: - mpd - the maximum dimension allowed for the projected problem
292: Options Database Keys:
293: + -qep_nev <nev> - Sets the number of eigenvalues
294: . -qep_ncv <ncv> - Sets the dimension of the subspace
295: - -qep_mpd <mpd> - Sets the maximum projected dimension
297: Notes:
298: Use PETSC_IGNORE to retain the previous value of any parameter.
300: Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
301: dependent on the solution method.
303: The parameters ncv and mpd are intimately related, so that the user is advised
304: to set one of them at most. Normal usage is the following
305: + - In cases where nev is small, the user sets ncv (a reasonable default is 2*nev).
306: - - In cases where nev is large, the user sets mpd.
308: The value of ncv should always be between nev and (nev+mpd), typically
309: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
310: a smaller value should be used.
312: Level: intermediate
314: .seealso: QEPGetDimensions()
315: @*/
316: PetscErrorCode QEPSetDimensions(QEP qep,PetscInt nev,PetscInt ncv,PetscInt mpd)
317: {
321: if( nev != PETSC_IGNORE ) {
322: if (nev<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
323: qep->nev = nev;
324: qep->setupcalled = 0;
325: }
326: if( ncv != PETSC_IGNORE ) {
327: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
328: qep->ncv = 0;
329: } else {
330: if (ncv<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
331: qep->ncv = ncv;
332: }
333: qep->setupcalled = 0;
334: }
335: if( mpd != PETSC_IGNORE ) {
336: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
337: qep->mpd = 0;
338: } else {
339: if (mpd<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
340: qep->mpd = mpd;
341: }
342: }
343: return(0);
344: }
348: /*@
349: QEPSetWhichEigenpairs - Specifies which portion of the spectrum is
350: to be sought.
352: Collective on QEP
354: Input Parameters:
355: + qep - eigensolver context obtained from QEPCreate()
356: - which - the portion of the spectrum to be sought
358: Possible values:
359: The parameter 'which' can have one of these values
360:
361: + QEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
362: . QEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
363: . QEP_LARGEST_REAL - largest real parts
364: . QEP_SMALLEST_REAL - smallest real parts
365: . QEP_LARGEST_IMAGINARY - largest imaginary parts
366: - QEP_SMALLEST_IMAGINARY - smallest imaginary parts
368: Options Database Keys:
369: + -qep_largest_magnitude - Sets largest eigenvalues in magnitude
370: . -qep_smallest_magnitude - Sets smallest eigenvalues in magnitude
371: . -qep_largest_real - Sets largest real parts
372: . -qep_smallest_real - Sets smallest real parts
373: . -qep_largest_imaginary - Sets largest imaginary parts
374: - -qep_smallest_imaginary - Sets smallest imaginary parts
376: Notes:
377: Not all eigensolvers implemented in QEP account for all the possible values
378: stated above. If SLEPc is compiled for real numbers QEP_LARGEST_IMAGINARY
379: and QEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
380: for eigenvalue selection.
381:
382: Level: intermediate
384: .seealso: QEPGetWhichEigenpairs(), QEPWhich
385: @*/
386: PetscErrorCode QEPSetWhichEigenpairs(QEP qep,QEPWhich which)
387: {
390: if (which!=PETSC_IGNORE) {
391: if (which==PETSC_DECIDE || which==PETSC_DEFAULT) qep->which = (QEPWhich)0;
392: else switch (which) {
393: case QEP_LARGEST_MAGNITUDE:
394: case QEP_SMALLEST_MAGNITUDE:
395: case QEP_LARGEST_REAL:
396: case QEP_SMALLEST_REAL:
397: case QEP_LARGEST_IMAGINARY:
398: case QEP_SMALLEST_IMAGINARY:
399: if (qep->which != which) {
400: qep->setupcalled = 0;
401: qep->which = which;
402: }
403: break;
404: default:
405: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
406: }
407: }
408: return(0);
409: }
413: /*@C
414: QEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
415: sought.
417: Not Collective
419: Input Parameter:
420: . qep - eigensolver context obtained from QEPCreate()
422: Output Parameter:
423: . which - the portion of the spectrum to be sought
425: Notes:
426: See QEPSetWhichEigenpairs() for possible values of 'which'.
428: Level: intermediate
430: .seealso: QEPSetWhichEigenpairs(), QEPWhich
431: @*/
432: PetscErrorCode QEPGetWhichEigenpairs(QEP qep,QEPWhich *which)
433: {
437: *which = qep->which;
438: return(0);
439: }
443: /*@
444: QEPSetLeftVectorsWanted - Specifies which eigenvectors are required.
446: Collective on QEP
448: Input Parameters:
449: + qep - the quadratic eigensolver context
450: - leftvecs - whether left eigenvectors are required or not
452: Options Database Keys:
453: . -qep_left_vectors <boolean> - Sets/resets the boolean flag 'leftvecs'
455: Notes:
456: If the user sets leftvecs=PETSC_TRUE then the solver uses a variant of
457: the algorithm that computes both right and left eigenvectors. This is
458: usually much more costly. This option is not available in all solvers.
460: Level: intermediate
462: .seealso: QEPGetLeftVectorsWanted(), QEPGetEigenvectorLeft()
463: @*/
464: PetscErrorCode QEPSetLeftVectorsWanted(QEP qep,PetscTruth leftvecs)
465: {
468: if (qep->leftvecs != leftvecs) {
469: qep->leftvecs = leftvecs;
470: qep->setupcalled = 0;
471: }
472: return(0);
473: }
477: /*@C
478: QEPGetLeftVectorsWanted - Returns the flag indicating whether left
479: eigenvectors are required or not.
481: Not Collective
483: Input Parameter:
484: . qep - the eigensolver context
486: Output Parameter:
487: . leftvecs - the returned flag
489: Level: intermediate
491: .seealso: QEPSetLeftVectorsWanted()
492: @*/
493: PetscErrorCode QEPGetLeftVectorsWanted(QEP qep,PetscTruth *leftvecs)
494: {
498: *leftvecs = qep->leftvecs;
499: return(0);
500: }
504: /*@
505: QEPGetScaleFactor - Gets the factor used for scaling the quadratic eigenproblem.
507: Not Collective
509: Input Parameter:
510: . qep - the quadratic eigensolver context
511:
512: Output Parameters:
513: . alpha - the scaling factor
515: Notes:
516: If the user did not specify a scaling factor, then after QEPSolve() the
517: default value is returned.
519: Level: intermediate
521: .seealso: QEPSetScaleFactor(), QEPSolve()
522: @*/
523: PetscErrorCode QEPGetScaleFactor(QEP qep,PetscReal *alpha)
524: {
527: if (alpha) *alpha = qep->sfactor;
528: return(0);
529: }
533: /*@
534: QEPSetScaleFactor - Sets the scaling factor to be used for scaling the
535: quadratic problem before attempting to solve.
537: Collective on QEP
539: Input Parameters:
540: + qep - the quadratic eigensolver context
541: - alpha - the scaling factor
543: Options Database Keys:
544: . -qep_scale <alpha> - Sets the scaling factor
546: Notes:
547: For the problem (l^2*M + l*C + K)*x = 0, the effect of scaling is to work
548: with matrices (alpha^2*M, alpha*C, K), then scale the computed eigenvalue.
550: The default is to scale with alpha = norm(K)/norm(M).
552: Level: intermediate
554: .seealso: QEPGetScaleFactor()
555: @*/
556: PetscErrorCode QEPSetScaleFactor(QEP qep,PetscReal alpha)
557: {
560: if (alpha != PETSC_IGNORE) {
561: if (alpha == PETSC_DEFAULT || alpha == PETSC_DECIDE) {
562: qep->sfactor = 0.0;
563: } else {
564: if (alpha < 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of alpha. Must be > 0");
565: qep->sfactor = alpha;
566: }
567: }
568: return(0);
569: }
573: /*@
574: QEPSetProblemType - Specifies the type of the quadratic eigenvalue problem.
576: Collective on QEP
578: Input Parameters:
579: + qep - the quadratic eigensolver context
580: - type - a known type of quadratic eigenvalue problem
582: Options Database Keys:
583: + -qep_general - general problem with no particular structure
584: . -qep_hermitian - problem whose coefficient matrices are Hermitian
585: - -qep_gyroscopic - problem with Hamiltonian structure
586:
587: Notes:
588: Allowed values for the problem type are: general (QEP_GENERAL), Hermitian
589: (QEP_HERMITIAN), and gyroscopic (QEP_GYROSCOPIC).
591: This function is used to instruct SLEPc to exploit certain structure in
592: the quadratic eigenproblem. By default, no particular structure is assumed.
594: If the problem matrices are Hermitian (symmetric in the real case) or
595: Hermitian/skew-Hermitian then the solver can exploit this fact to perform
596: less operations or provide better stability.
598: Level: intermediate
600: .seealso: QEPSetOperators(), QEPSetType(), QEPGetProblemType(), QEPProblemType
601: @*/
602: PetscErrorCode QEPSetProblemType(QEP qep,QEPProblemType type)
603: {
606: if (type!=QEP_GENERAL && type!=QEP_HERMITIAN && type!=QEP_GYROSCOPIC)
607: SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
608: qep->problem_type = type;
609: return(0);
610: }
614: /*@C
615: QEPGetProblemType - Gets the problem type from the QEP object.
617: Not Collective
619: Input Parameter:
620: . qep - the quadratic eigensolver context
622: Output Parameter:
623: . type - name of QEP problem type
625: Level: intermediate
627: .seealso: QEPSetProblemType(), QEPProblemType
628: @*/
629: PetscErrorCode QEPGetProblemType(QEP qep,QEPProblemType *type)
630: {
634: *type = qep->problem_type;
635: return(0);
636: }
640: /*@C
641: QEPSetConvergenceTest - Sets a function to compute the error estimate used in
642: the convergence test.
644: Collective on QEP
646: Input Parameters:
647: + qep - eigensolver context obtained from QEPCreate()
648: . func - a pointer to the convergence test function
649: - ctx - a context pointer (the last parameter to the convergence test function)
651: Calling Sequence of func:
652: $ func(QEP qep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal* errest,void *ctx)
654: + qep - eigensolver context obtained from QEPCreate()
655: . eigr - real part of the eigenvalue
656: . eigi - imaginary part of the eigenvalue
657: . res - residual norm associated to the eigenpair
658: . errest - (output) computed error estimate
659: - ctx - optional context, as set by QEPSetConvergenceTest()
661: Note:
662: If the error estimate returned by the convergence test function is less than
663: the tolerance, then the eigenvalue is accepted as converged.
665: Level: advanced
667: .seealso: QEPSetTolerances()
668: @*/
669: EXTERN PetscErrorCode QEPSetConvergenceTest(QEP qep,PetscErrorCode (*func)(QEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx)
670: {
672: qep->conv_func = func;
673: qep->conv_ctx = ctx;
674: return(0);
675: }
679: /*@
680: QEPSetTrackAll - Specifies if the solver must compute the residual of all
681: approximate eigenpairs or not.
683: Collective on QEP
685: Input Parameters:
686: + qep - the eigensolver context
687: - trackall - whether compute all residuals or not
689: Notes:
690: If the user sets trackall=PETSC_TRUE then the solver explicitly computes
691: the residual for each eigenpair approximation. Computing the residual is
692: usually an expensive operation and solvers commonly compute the associated
693: residual to the first unconverged eigenpair.
694: The options '-qep_monitor_all' and '-qep_monitor_draw_all' automatically
695: activates this option.
697: Level: intermediate
699: .seealso: EPSGetTrackAll()
700: @*/
701: PetscErrorCode QEPSetTrackAll(QEP qep,PetscTruth trackall)
702: {
705: qep->trackall = trackall;
706: return(0);
707: }
711: /*@
712: QEPGetTrackAll - Returns the flag indicating whether all residual norms must
713: be computed or not.
715: Not Collective
717: Input Parameter:
718: . qep - the eigensolver context
720: Output Parameter:
721: . trackall - the returned flag
723: Level: intermediate
725: .seealso: EPSSetTrackAll()
726: @*/
727: PetscErrorCode QEPGetTrackAll(QEP qep,PetscTruth *trackall)
728: {
732: *trackall = qep->trackall;
733: return(0);
734: }
738: /*@C
739: QEPSetOptionsPrefix - Sets the prefix used for searching for all
740: QEP options in the database.
742: Collective on QEP
744: Input Parameters:
745: + qep - the quadratic eigensolver context
746: - prefix - the prefix string to prepend to all QEP option requests
748: Notes:
749: A hyphen (-) must NOT be given at the beginning of the prefix name.
750: The first character of all runtime options is AUTOMATICALLY the
751: hyphen.
753: For example, to distinguish between the runtime options for two
754: different QEP contexts, one could call
755: .vb
756: QEPSetOptionsPrefix(qep1,"qeig1_")
757: QEPSetOptionsPrefix(qep2,"qeig2_")
758: .ve
760: Level: advanced
762: .seealso: QEPAppendOptionsPrefix(), QEPGetOptionsPrefix()
763: @*/
764: PetscErrorCode QEPSetOptionsPrefix(QEP qep,const char *prefix)
765: {
769: PetscObjectSetOptionsPrefix((PetscObject)qep,prefix);
770: IPSetOptionsPrefix(qep->ip,prefix);
771: IPAppendOptionsPrefix(qep->ip,"qep_");
772: return(0);
773: }
774:
777: /*@C
778: QEPAppendOptionsPrefix - Appends to the prefix used for searching for all
779: QEP options in the database.
781: Collective on QEP
783: Input Parameters:
784: + qep - the quadratic eigensolver context
785: - prefix - the prefix string to prepend to all QEP option requests
787: Notes:
788: A hyphen (-) must NOT be given at the beginning of the prefix name.
789: The first character of all runtime options is AUTOMATICALLY the hyphen.
791: Level: advanced
793: .seealso: QEPSetOptionsPrefix(), QEPGetOptionsPrefix()
794: @*/
795: PetscErrorCode QEPAppendOptionsPrefix(QEP qep,const char *prefix)
796: {
800: PetscObjectAppendOptionsPrefix((PetscObject)qep, prefix);
801: IPSetOptionsPrefix(qep->ip,prefix);
802: IPAppendOptionsPrefix(qep->ip,"qep_");
803: return(0);
804: }
808: /*@C
809: QEPGetOptionsPrefix - Gets the prefix used for searching for all
810: QEP options in the database.
812: Not Collective
814: Input Parameters:
815: . qep - the quadratic eigensolver context
817: Output Parameters:
818: . prefix - pointer to the prefix string used is returned
820: Notes: On the fortran side, the user should pass in a string 'prefix' of
821: sufficient length to hold the prefix.
823: Level: advanced
825: .seealso: QEPSetOptionsPrefix(), QEPAppendOptionsPrefix()
826: @*/
827: PetscErrorCode QEPGetOptionsPrefix(QEP qep,const char *prefix[])
828: {
833: PetscObjectGetOptionsPrefix((PetscObject)qep,prefix);
834: return(0);
835: }