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-2013, 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/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: PetscScalar s;
51: PetscInt i,j,k;
52: PetscViewer monviewer;
53: SlepcConvMonitor ctx;
57: if (!QEPRegisterAllCalled) { QEPRegisterAll(); }
58: PetscObjectOptionsBegin((PetscObject)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 = 0;
74: PetscOptionsReal("-qep_scale","Scale factor","QEPSetScaleFactor",qep->sfactor,&r,NULL);
75: QEPSetScaleFactor(qep,r);
77: r = i = 0;
78: PetscOptionsInt("-qep_max_it","Maximum number of iterations","QEPSetTolerances",qep->max_it,&i,NULL);
79: PetscOptionsReal("-qep_tol","Tolerance","QEPSetTolerances",qep->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:qep->tol,&r,NULL);
80: QEPSetTolerances(qep,r,i);
81: PetscOptionsBoolGroupBegin("-qep_convergence_default","Default (relative error) convergence test","QEPSetConvergenceTest",&flg);
82: if (flg) { QEPSetConvergenceTest(qep,QEPConvergedDefault,NULL); }
83: PetscOptionsBoolGroupEnd("-qep_convergence_absolute","Absolute error convergence test","QEPSetConvergenceTest",&flg);
84: if (flg) { QEPSetConvergenceTest(qep,QEPConvergedAbsolute,NULL); }
86: i = j = k = 0;
87: PetscOptionsInt("-qep_nev","Number of eigenvalues to compute","QEPSetDimensions",qep->nev,&i,NULL);
88: PetscOptionsInt("-qep_ncv","Number of basis vectors","QEPSetDimensions",qep->ncv,&j,NULL);
89: PetscOptionsInt("-qep_mpd","Maximum dimension of projected problem","QEPSetDimensions",qep->mpd,&k,NULL);
90: QEPSetDimensions(qep,i,j,k);
92: PetscOptionsScalar("-qep_target","Value of the target","QEPSetTarget",qep->target,&s,&flg);
93: if (flg) {
94: QEPSetWhichEigenpairs(qep,QEP_TARGET_MAGNITUDE);
95: QEPSetTarget(qep,s);
96: }
98: /* -----------------------------------------------------------------------*/
99: /*
100: Cancels all monitors hardwired into code before call to QEPSetFromOptions()
101: */
102: flg = PETSC_FALSE;
103: PetscOptionsBool("-qep_monitor_cancel","Remove any hardwired monitor routines","QEPMonitorCancel",flg,&flg,NULL);
104: if (flg) {
105: QEPMonitorCancel(qep);
106: }
107: /*
108: Prints approximate eigenvalues and error estimates at each iteration
109: */
110: PetscOptionsString("-qep_monitor","Monitor first unconverged approximate eigenvalue and error estimate","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
111: if (flg) {
112: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)qep),monfilename,&monviewer);
113: QEPMonitorSet(qep,QEPMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
114: }
115: PetscOptionsString("-qep_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
116: if (flg) {
117: PetscNew(struct _n_SlepcConvMonitor,&ctx);
118: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)qep),monfilename,&ctx->viewer);
119: QEPMonitorSet(qep,QEPMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
120: }
121: PetscOptionsString("-qep_monitor_all","Monitor approximate eigenvalues and error estimates","QEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
122: if (flg) {
123: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)qep),monfilename,&monviewer);
124: QEPMonitorSet(qep,QEPMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
125: QEPSetTrackAll(qep,PETSC_TRUE);
126: }
127: flg = PETSC_FALSE;
128: PetscOptionsBool("-qep_monitor_lg","Monitor first unconverged approximate error estimate graphically","QEPMonitorSet",flg,&flg,NULL);
129: if (flg) {
130: QEPMonitorSet(qep,QEPMonitorLG,NULL,NULL);
131: }
132: flg = PETSC_FALSE;
133: PetscOptionsBool("-qep_monitor_lg_all","Monitor error estimates graphically","QEPMonitorSet",flg,&flg,NULL);
134: if (flg) {
135: QEPMonitorSet(qep,QEPMonitorLGAll,NULL,NULL);
136: QEPSetTrackAll(qep,PETSC_TRUE);
137: }
138: /* -----------------------------------------------------------------------*/
140: PetscOptionsBoolGroupBegin("-qep_largest_magnitude","compute largest eigenvalues in magnitude","QEPSetWhichEigenpairs",&flg);
141: if (flg) { QEPSetWhichEigenpairs(qep,QEP_LARGEST_MAGNITUDE); }
142: PetscOptionsBoolGroup("-qep_smallest_magnitude","compute smallest eigenvalues in magnitude","QEPSetWhichEigenpairs",&flg);
143: if (flg) { QEPSetWhichEigenpairs(qep,QEP_SMALLEST_MAGNITUDE); }
144: PetscOptionsBoolGroup("-qep_largest_real","compute largest real parts","QEPSetWhichEigenpairs",&flg);
145: if (flg) { QEPSetWhichEigenpairs(qep,QEP_LARGEST_REAL); }
146: PetscOptionsBoolGroup("-qep_smallest_real","compute smallest real parts","QEPSetWhichEigenpairs",&flg);
147: if (flg) { QEPSetWhichEigenpairs(qep,QEP_SMALLEST_REAL); }
148: PetscOptionsBoolGroup("-qep_largest_imaginary","compute largest imaginary parts","QEPSetWhichEigenpairs",&flg);
149: if (flg) { QEPSetWhichEigenpairs(qep,QEP_LARGEST_IMAGINARY); }
150: PetscOptionsBoolGroup("-qep_smallest_imaginary","compute smallest imaginary parts","QEPSetWhichEigenpairs",&flg);
151: if (flg) { QEPSetWhichEigenpairs(qep,QEP_SMALLEST_IMAGINARY); }
152: PetscOptionsBoolGroup("-qep_target_magnitude","compute nearest eigenvalues to target","QEPSetWhichEigenpairs",&flg);
153: if (flg) { QEPSetWhichEigenpairs(qep,QEP_TARGET_MAGNITUDE); }
154: PetscOptionsBoolGroup("-qep_target_real","compute eigenvalues with real parts close to target","QEPSetWhichEigenpairs",&flg);
155: if (flg) { QEPSetWhichEigenpairs(qep,QEP_TARGET_REAL); }
156: PetscOptionsBoolGroupEnd("-qep_target_imaginary","compute eigenvalues with imaginary parts close to target","QEPSetWhichEigenpairs",&flg);
157: if (flg) { QEPSetWhichEigenpairs(qep,QEP_TARGET_IMAGINARY); }
159: PetscOptionsBool("-qep_left_vectors","Compute left eigenvectors also","QEPSetLeftVectorsWanted",qep->leftvecs,&val,&flg);
160: if (flg) {
161: QEPSetLeftVectorsWanted(qep,val);
162: }
164: PetscOptionsName("-qep_view","Print detailed information on solver used","QEPView",0);
165: PetscOptionsName("-qep_plot_eigs","Make a plot of the computed eigenvalues","QEPSolve",0);
167: if (qep->ops->setfromoptions) {
168: (*qep->ops->setfromoptions)(qep);
169: }
170: PetscObjectProcessOptionsHandlers((PetscObject)qep);
171: PetscOptionsEnd();
173: if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
174: IPSetFromOptions(qep->ip);
175: if (!qep->ds) { QEPGetDS(qep,&qep->ds); }
176: DSSetFromOptions(qep->ds);
177: if (!qep->st) { QEPGetST(qep,&qep->st); }
178: STSetFromOptions(qep->st);
179: PetscRandomSetFromOptions(qep->rand);
180: return(0);
181: }
185: /*@
186: QEPGetTolerances - Gets the tolerance and maximum iteration count used
187: by the QEP convergence tests.
189: Not Collective
191: Input Parameter:
192: . qep - the quadratic eigensolver context
194: Output Parameters:
195: + tol - the convergence tolerance
196: - maxits - maximum number of iterations
198: Notes:
199: The user can specify NULL for any parameter that is not needed.
201: Level: intermediate
203: .seealso: QEPSetTolerances()
204: @*/
205: PetscErrorCode QEPGetTolerances(QEP qep,PetscReal *tol,PetscInt *maxits)
206: {
209: if (tol) *tol = qep->tol;
210: if (maxits) *maxits = qep->max_it;
211: return(0);
212: }
216: /*@
217: QEPSetTolerances - Sets the tolerance and maximum iteration count used
218: by the QEP convergence tests.
220: Logically Collective on QEP
222: Input Parameters:
223: + qep - the quadratic eigensolver context
224: . tol - the convergence tolerance
225: - maxits - maximum number of iterations to use
227: Options Database Keys:
228: + -qep_tol <tol> - Sets the convergence tolerance
229: - -qep_max_it <maxits> - Sets the maximum number of iterations allowed
231: Notes:
232: Pass 0 for an argument that need not be changed.
234: Use PETSC_DECIDE for maxits to assign a reasonably good value, which is
235: dependent on the solution method.
237: Level: intermediate
239: .seealso: QEPGetTolerances()
240: @*/
241: PetscErrorCode QEPSetTolerances(QEP qep,PetscReal tol,PetscInt maxits)
242: {
247: if (tol) {
248: if (tol == PETSC_DEFAULT) {
249: qep->tol = PETSC_DEFAULT;
250: } else {
251: if (tol < 0.0) SETERRQ(PetscObjectComm((PetscObject)qep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
252: qep->tol = tol;
253: }
254: }
255: if (maxits) {
256: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
257: qep->max_it = 0;
258: qep->setupcalled = 0;
259: } else {
260: if (maxits < 0) SETERRQ(PetscObjectComm((PetscObject)qep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
261: qep->max_it = maxits;
262: }
263: }
264: return(0);
265: }
269: /*@
270: QEPGetDimensions - Gets the number of eigenvalues to compute
271: and the dimension of the subspace.
273: Not Collective
275: Input Parameter:
276: . qep - the quadratic eigensolver context
278: Output Parameters:
279: + nev - number of eigenvalues to compute
280: . ncv - the maximum dimension of the subspace to be used by the solver
281: - mpd - the maximum dimension allowed for the projected problem
283: Notes:
284: The user can specify NULL for any parameter that is not needed.
286: Level: intermediate
288: .seealso: QEPSetDimensions()
289: @*/
290: PetscErrorCode QEPGetDimensions(QEP qep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
291: {
294: if (nev) *nev = qep->nev;
295: if (ncv) *ncv = qep->ncv;
296: if (mpd) *mpd = qep->mpd;
297: return(0);
298: }
302: /*@
303: QEPSetDimensions - Sets the number of eigenvalues to compute
304: and the dimension of the subspace.
306: Logically Collective on QEP
308: Input Parameters:
309: + qep - the quadratic eigensolver context
310: . nev - number of eigenvalues to compute
311: . ncv - the maximum dimension of the subspace to be used by the solver
312: - mpd - the maximum dimension allowed for the projected problem
314: Options Database Keys:
315: + -qep_nev <nev> - Sets the number of eigenvalues
316: . -qep_ncv <ncv> - Sets the dimension of the subspace
317: - -qep_mpd <mpd> - Sets the maximum projected dimension
319: Notes:
320: Pass 0 to retain the previous value of any parameter.
322: Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
323: dependent on the solution method.
325: The parameters ncv and mpd are intimately related, so that the user is advised
326: to set one of them at most. Normal usage is that
327: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
328: (b) in cases where nev is large, the user sets mpd.
330: The value of ncv should always be between nev and (nev+mpd), typically
331: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
332: a smaller value should be used.
334: Level: intermediate
336: .seealso: QEPGetDimensions()
337: @*/
338: PetscErrorCode QEPSetDimensions(QEP qep,PetscInt nev,PetscInt ncv,PetscInt mpd)
339: {
345: if (nev) {
346: if (nev<1) SETERRQ(PetscObjectComm((PetscObject)qep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
347: qep->nev = nev;
348: qep->setupcalled = 0;
349: }
350: if (ncv) {
351: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
352: qep->ncv = 0;
353: } else {
354: if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)qep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
355: qep->ncv = ncv;
356: }
357: qep->setupcalled = 0;
358: }
359: if (mpd) {
360: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
361: qep->mpd = 0;
362: } else {
363: if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)qep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
364: qep->mpd = mpd;
365: }
366: }
367: return(0);
368: }
372: /*@
373: QEPSetWhichEigenpairs - Specifies which portion of the spectrum is
374: to be sought.
376: Logically Collective on QEP
378: Input Parameters:
379: + qep - eigensolver context obtained from QEPCreate()
380: - which - the portion of the spectrum to be sought
382: Possible values:
383: The parameter 'which' can have one of these values
385: + QEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
386: . QEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
387: . QEP_LARGEST_REAL - largest real parts
388: . QEP_SMALLEST_REAL - smallest real parts
389: . QEP_LARGEST_IMAGINARY - largest imaginary parts
390: . QEP_SMALLEST_IMAGINARY - smallest imaginary parts
391: . QEP_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
392: . QEP_TARGET_REAL - eigenvalues with real part closest to target
393: - QEP_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
395: Options Database Keys:
396: + -qep_largest_magnitude - Sets largest eigenvalues in magnitude
397: . -qep_smallest_magnitude - Sets smallest eigenvalues in magnitude
398: . -qep_largest_real - Sets largest real parts
399: . -qep_smallest_real - Sets smallest real parts
400: . -qep_largest_imaginary - Sets largest imaginary parts
401: . -qep_smallest_imaginary - Sets smallest imaginary parts
402: . -qep_target_magnitude - Sets eigenvalues closest to target
403: . -qep_target_real - Sets real parts closest to target
404: - -qep_target_imaginary - Sets imaginary parts closest to target
406: Notes:
407: Not all eigensolvers implemented in QEP account for all the possible values
408: stated above. If SLEPc is compiled for real numbers QEP_LARGEST_IMAGINARY
409: and QEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
410: for eigenvalue selection.
412: Level: intermediate
414: .seealso: QEPGetWhichEigenpairs(), QEPWhich
415: @*/
416: PetscErrorCode QEPSetWhichEigenpairs(QEP qep,QEPWhich which)
417: {
421: if (which) {
422: if (which==PETSC_DECIDE || which==PETSC_DEFAULT) qep->which = (QEPWhich)0;
423: else switch (which) {
424: case QEP_LARGEST_MAGNITUDE:
425: case QEP_SMALLEST_MAGNITUDE:
426: case QEP_LARGEST_REAL:
427: case QEP_SMALLEST_REAL:
428: case QEP_LARGEST_IMAGINARY:
429: case QEP_SMALLEST_IMAGINARY:
430: case QEP_TARGET_MAGNITUDE:
431: case QEP_TARGET_REAL:
432: #if defined(PETSC_USE_COMPLEX)
433: case QEP_TARGET_IMAGINARY:
434: #endif
435: if (qep->which != which) {
436: qep->setupcalled = 0;
437: qep->which = which;
438: }
439: break;
440: default:
441: SETERRQ(PetscObjectComm((PetscObject)qep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
442: }
443: }
444: return(0);
445: }
449: /*@C
450: QEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
451: sought.
453: Not Collective
455: Input Parameter:
456: . qep - eigensolver context obtained from QEPCreate()
458: Output Parameter:
459: . which - the portion of the spectrum to be sought
461: Notes:
462: See QEPSetWhichEigenpairs() for possible values of 'which'.
464: Level: intermediate
466: .seealso: QEPSetWhichEigenpairs(), QEPWhich
467: @*/
468: PetscErrorCode QEPGetWhichEigenpairs(QEP qep,QEPWhich *which)
469: {
473: *which = qep->which;
474: return(0);
475: }
479: /*@
480: QEPSetLeftVectorsWanted - Specifies which eigenvectors are required.
482: Logically Collective on QEP
484: Input Parameters:
485: + qep - the quadratic eigensolver context
486: - leftvecs - whether left eigenvectors are required or not
488: Options Database Keys:
489: . -qep_left_vectors <boolean> - Sets/resets the boolean flag 'leftvecs'
491: Notes:
492: If the user sets leftvecs=PETSC_TRUE then the solver uses a variant of
493: the algorithm that computes both right and left eigenvectors. This is
494: usually much more costly. This option is not available in all solvers.
496: Level: intermediate
498: .seealso: QEPGetLeftVectorsWanted()
499: @*/
500: PetscErrorCode QEPSetLeftVectorsWanted(QEP qep,PetscBool leftvecs)
501: {
505: if (qep->leftvecs != leftvecs) {
506: qep->leftvecs = leftvecs;
507: qep->setupcalled = 0;
508: }
509: return(0);
510: }
514: /*@C
515: QEPGetLeftVectorsWanted - Returns the flag indicating whether left
516: eigenvectors are required or not.
518: Not Collective
520: Input Parameter:
521: . qep - the eigensolver context
523: Output Parameter:
524: . leftvecs - the returned flag
526: Level: intermediate
528: .seealso: QEPSetLeftVectorsWanted()
529: @*/
530: PetscErrorCode QEPGetLeftVectorsWanted(QEP qep,PetscBool *leftvecs)
531: {
535: *leftvecs = qep->leftvecs;
536: return(0);
537: }
541: /*@
542: QEPGetScaleFactor - Gets the factor used for scaling the quadratic eigenproblem.
544: Not Collective
546: Input Parameter:
547: . qep - the quadratic eigensolver context
549: Output Parameters:
550: . alpha - the scaling factor
552: Notes:
553: If the user did not specify a scaling factor, then after QEPSolve() the
554: default value is returned.
556: Level: intermediate
558: .seealso: QEPSetScaleFactor(), QEPSolve()
559: @*/
560: PetscErrorCode QEPGetScaleFactor(QEP qep,PetscReal *alpha)
561: {
565: *alpha = qep->sfactor;
566: return(0);
567: }
571: /*@
572: QEPSetScaleFactor - Sets the scaling factor to be used for scaling the
573: quadratic problem before attempting to solve.
575: Logically Collective on QEP
577: Input Parameters:
578: + qep - the quadratic eigensolver context
579: - alpha - the scaling factor
581: Options Database Keys:
582: . -qep_scale <alpha> - Sets the scaling factor
584: Notes:
585: For the problem (l^2*M + l*C + K)*x = 0, the effect of scaling is to work
586: with matrices (alpha^2*M, alpha*C, K), then scale the computed eigenvalue.
588: The default is to scale with alpha = norm(K)/norm(M).
590: Level: intermediate
592: .seealso: QEPGetScaleFactor()
593: @*/
594: PetscErrorCode QEPSetScaleFactor(QEP qep,PetscReal alpha)
595: {
599: if (alpha) {
600: if (alpha == PETSC_DEFAULT || alpha == PETSC_DECIDE) {
601: qep->sfactor = 0.0;
602: qep->sfactor_set = PETSC_FALSE;
603: } else {
604: if (alpha < 0.0) SETERRQ(PetscObjectComm((PetscObject)qep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of alpha. Must be > 0");
605: qep->sfactor = alpha;
606: qep->sfactor_set = PETSC_TRUE;
607: }
608: }
609: return(0);
610: }
614: /*@
615: QEPSetProblemType - Specifies the type of the quadratic eigenvalue problem.
617: Logically Collective on QEP
619: Input Parameters:
620: + qep - the quadratic eigensolver context
621: - type - a known type of quadratic eigenvalue problem
623: Options Database Keys:
624: + -qep_general - general problem with no particular structure
625: . -qep_hermitian - problem whose coefficient matrices are Hermitian
626: - -qep_gyroscopic - problem with Hamiltonian structure
628: Notes:
629: Allowed values for the problem type are: general (QEP_GENERAL), Hermitian
630: (QEP_HERMITIAN), and gyroscopic (QEP_GYROSCOPIC).
632: This function is used to instruct SLEPc to exploit certain structure in
633: the quadratic eigenproblem. By default, no particular structure is assumed.
635: If the problem matrices are Hermitian (symmetric in the real case) or
636: Hermitian/skew-Hermitian then the solver can exploit this fact to perform
637: less operations or provide better stability.
639: Level: intermediate
641: .seealso: QEPSetOperators(), QEPSetType(), QEPGetProblemType(), QEPProblemType
642: @*/
643: PetscErrorCode QEPSetProblemType(QEP qep,QEPProblemType type)
644: {
648: if (type!=QEP_GENERAL && type!=QEP_HERMITIAN && type!=QEP_GYROSCOPIC)
649: SETERRQ(PetscObjectComm((PetscObject)qep),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
650: qep->problem_type = type;
651: return(0);
652: }
656: /*@C
657: QEPGetProblemType - Gets the problem type from the QEP object.
659: Not Collective
661: Input Parameter:
662: . qep - the quadratic eigensolver context
664: Output Parameter:
665: . type - name of QEP problem type
667: Level: intermediate
669: .seealso: QEPSetProblemType(), QEPProblemType
670: @*/
671: PetscErrorCode QEPGetProblemType(QEP qep,QEPProblemType *type)
672: {
676: *type = qep->problem_type;
677: return(0);
678: }
682: /*@C
683: QEPSetConvergenceTest - Sets a function to compute the error estimate used in
684: the convergence test.
686: Logically Collective on QEP
688: Input Parameters:
689: + qep - eigensolver context obtained from QEPCreate()
690: . func - a pointer to the convergence test function
691: - ctx - a context pointer (the last parameter to the convergence test function)
693: Calling Sequence of func:
694: $ func(QEP qep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal* errest,void *ctx)
696: + qep - eigensolver context obtained from QEPCreate()
697: . eigr - real part of the eigenvalue
698: . eigi - imaginary part of the eigenvalue
699: . res - residual norm associated to the eigenpair
700: . errest - (output) computed error estimate
701: - ctx - optional context, as set by QEPSetConvergenceTest()
703: Note:
704: If the error estimate returned by the convergence test function is less than
705: the tolerance, then the eigenvalue is accepted as converged.
707: Level: advanced
709: .seealso: QEPSetTolerances()
710: @*/
711: PetscErrorCode QEPSetConvergenceTest(QEP qep,PetscErrorCode (*func)(QEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx)
712: {
715: qep->converged = func;
716: qep->convergedctx = ctx;
717: return(0);
718: }
722: /*@
723: QEPSetTrackAll - Specifies if the solver must compute the residual of all
724: approximate eigenpairs or not.
726: Logically Collective on QEP
728: Input Parameters:
729: + qep - the eigensolver context
730: - trackall - whether compute all residuals or not
732: Notes:
733: If the user sets trackall=PETSC_TRUE then the solver explicitly computes
734: the residual for each eigenpair approximation. Computing the residual is
735: usually an expensive operation and solvers commonly compute the associated
736: residual to the first unconverged eigenpair.
737: The options '-qep_monitor_all' and '-qep_monitor_lg_all' automatically
738: activate this option.
740: Level: intermediate
742: .seealso: QEPGetTrackAll()
743: @*/
744: PetscErrorCode QEPSetTrackAll(QEP qep,PetscBool trackall)
745: {
749: qep->trackall = trackall;
750: return(0);
751: }
755: /*@
756: QEPGetTrackAll - Returns the flag indicating whether all residual norms must
757: be computed or not.
759: Not Collective
761: Input Parameter:
762: . qep - the eigensolver context
764: Output Parameter:
765: . trackall - the returned flag
767: Level: intermediate
769: .seealso: QEPSetTrackAll()
770: @*/
771: PetscErrorCode QEPGetTrackAll(QEP qep,PetscBool *trackall)
772: {
776: *trackall = qep->trackall;
777: return(0);
778: }
782: /*@C
783: QEPSetOptionsPrefix - Sets the prefix used for searching for all
784: QEP options in the database.
786: Logically Collective on QEP
788: Input Parameters:
789: + qep - the quadratic eigensolver context
790: - prefix - the prefix string to prepend to all QEP option requests
792: Notes:
793: A hyphen (-) must NOT be given at the beginning of the prefix name.
794: The first character of all runtime options is AUTOMATICALLY the
795: hyphen.
797: For example, to distinguish between the runtime options for two
798: different QEP contexts, one could call
799: .vb
800: QEPSetOptionsPrefix(qep1,"qeig1_")
801: QEPSetOptionsPrefix(qep2,"qeig2_")
802: .ve
804: Level: advanced
806: .seealso: QEPAppendOptionsPrefix(), QEPGetOptionsPrefix()
807: @*/
808: PetscErrorCode QEPSetOptionsPrefix(QEP qep,const char *prefix)
809: {
814: if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
815: IPSetOptionsPrefix(qep->ip,prefix);
816: if (!qep->ds) { QEPGetDS(qep,&qep->ds); }
817: DSSetOptionsPrefix(qep->ds,prefix);
818: PetscObjectSetOptionsPrefix((PetscObject)qep,prefix);
819: return(0);
820: }
824: /*@C
825: QEPAppendOptionsPrefix - Appends to the prefix used for searching for all
826: QEP options in the database.
828: Logically Collective on QEP
830: Input Parameters:
831: + qep - the quadratic eigensolver context
832: - prefix - the prefix string to prepend to all QEP option requests
834: Notes:
835: A hyphen (-) must NOT be given at the beginning of the prefix name.
836: The first character of all runtime options is AUTOMATICALLY the hyphen.
838: Level: advanced
840: .seealso: QEPSetOptionsPrefix(), QEPGetOptionsPrefix()
841: @*/
842: PetscErrorCode QEPAppendOptionsPrefix(QEP qep,const char *prefix)
843: {
845: PetscBool flg;
846: EPS eps;
850: if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
851: IPSetOptionsPrefix(qep->ip,prefix);
852: if (!qep->ds) { QEPGetDS(qep,&qep->ds); }
853: DSSetOptionsPrefix(qep->ds,prefix);
854: PetscObjectAppendOptionsPrefix((PetscObject)qep,prefix);
855: PetscObjectTypeCompare((PetscObject)qep,QEPLINEAR,&flg);
856: if (flg) {
857: QEPLinearGetEPS(qep,&eps);
858: EPSSetOptionsPrefix(eps,((PetscObject)qep)->prefix);
859: EPSAppendOptionsPrefix(eps,"qep_");
860: }
861: return(0);
862: }
866: /*@C
867: QEPGetOptionsPrefix - Gets the prefix used for searching for all
868: QEP options in the database.
870: Not Collective
872: Input Parameters:
873: . qep - the quadratic eigensolver context
875: Output Parameters:
876: . prefix - pointer to the prefix string used is returned
878: Notes: On the fortran side, the user should pass in a string 'prefix' of
879: sufficient length to hold the prefix.
881: Level: advanced
883: .seealso: QEPSetOptionsPrefix(), QEPAppendOptionsPrefix()
884: @*/
885: PetscErrorCode QEPGetOptionsPrefix(QEP qep,const char *prefix[])
886: {
892: PetscObjectGetOptionsPrefix((PetscObject)qep,prefix);
893: return(0);
894: }