Actual source code: opts.c
1: /*
2: EPS routines related to options that can be set via the command-line
3: or procedurally.
4: */
5: #include src/eps/epsimpl.h
9: /*@
10: EPSSetFromOptions - Sets EPS options from the options database.
11: This routine must be called before EPSSetUp() if the user is to be
12: allowed to set the solver type.
14: Collective on EPS
16: Input Parameters:
17: . eps - the eigensolver context
19: Notes:
20: To see all options, run your program with the -help option.
22: Level: beginner
24: .seealso:
25: @*/
26: PetscErrorCode EPSSetFromOptions(EPS eps)
27: {
29: char type[256],monfilename[PETSC_MAX_PATH_LEN];
30: PetscTruth flg;
31: const char *orth_list[3] = { "mgs" , "cgs", "ncgs" };
32: const char *ref_list[3] = { "never" , "ifneeded", "always" };
33: PetscReal eta;
34: PetscInt i,orth_type,ref_type;
35: PetscViewer monviewer;
39: PetscOptionsBegin(eps->comm,eps->prefix,"Eigenproblem Solver (EPS) Options","EPS");
40: PetscOptionsList("-eps_type","Eigenproblem Solver method","EPSSetType",EPSList,(char*)(eps->type_name?eps->type_name:EPSKRYLOVSCHUR),type,256,&flg);
41: if (flg) {
42: EPSSetType(eps,type);
43: }
45: PetscOptionsTruthGroupBegin("-eps_hermitian","hermitian eigenvalue problem","EPSSetProblemType",&flg);
46: if (flg) {EPSSetProblemType(eps,EPS_HEP);}
47: PetscOptionsTruthGroup("-eps_gen_hermitian","generalized hermitian eigenvalue problem","EPSSetProblemType",&flg);
48: if (flg) {EPSSetProblemType(eps,EPS_GHEP);}
49: PetscOptionsTruthGroup("-eps_non_hermitian","non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
50: if (flg) {EPSSetProblemType(eps,EPS_NHEP);}
51: PetscOptionsTruthGroupEnd("-eps_gen_non_hermitian","generalized non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
52: if (flg) {EPSSetProblemType(eps,EPS_GNHEP);}
54: /*
55: Set the type if it was never set.
56: */
57: if (!eps->type_name) {
58: EPSSetType(eps,EPSKRYLOVSCHUR);
59: }
61: PetscOptionsTruthGroupBegin("-eps_oneside","one-sided eigensolver","EPSSetClass",&flg);
62: if (flg) {EPSSetClass(eps,EPS_ONE_SIDE);}
63: PetscOptionsTruthGroupEnd("-eps_twoside","two-sided eigensolver","EPSSetClass",&flg);
64: if (flg) {EPSSetClass(eps,EPS_TWO_SIDE);}
66: orth_type = eps->orthog_type;
67: PetscOptionsEList("-eps_orthog_type","Orthogonalization method","EPSSetOrthogonalization",orth_list,3,orth_list[eps->orthog_type],&orth_type,&flg);
68: ref_type = eps->orthog_ref;
69: PetscOptionsEList("-eps_orthog_refinement","Iterative refinement mode during orthogonalization","EPSSetOrthogonalizationRefinement",ref_list,3,ref_list[eps->orthog_ref],&ref_type,&flg);
70: eta = eps->orthog_eta;
71: PetscOptionsReal("-eps_orthog_eta","Parameter of iterative refinement during orthogonalization","EPSSetOrthogonalizationRefinement",eta,&eta,PETSC_NULL);
72: EPSSetOrthogonalization(eps,(EPSOrthogonalizationType)orth_type,(EPSOrthogonalizationRefinementType)ref_type,eta);
74: PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&i,&flg);
75: if (flg) eps->max_it = i;
76: PetscOptionsReal("-eps_tol","Tolerance","KSPSetTolerances",eps->tol,&eps->tol,PETSC_NULL);
77: PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&i,&flg);
78: if (flg) {
79: if(i<1) SETERRQ(1,"Illegal value for option -eps_nev. Must be > 0");
80: eps->nev = i;
81: }
82: PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&i,&flg);
83: if (flg) {
84: if (i<1) SETERRQ(1,"Illegal value for option -eps_ncv. Must be > 0");
85: eps->ncv = i;
86: }
88: /* -----------------------------------------------------------------------*/
89: /*
90: Cancels all monitors hardwired into code before call to EPSSetFromOptions()
91: */
92: PetscOptionsName("-eps_cancelmonitors","Remove any hardwired monitor routines","EPSClearMonitor",&flg);
93: if (flg) {
94: EPSClearMonitor(eps);
95: }
96: /*
97: Prints approximate eigenvalues and error estimates at each iteration
98: */
99: PetscOptionsString("-eps_monitor","Monitor approximate eigenvalues and error estimates","EPSSetMonitor","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
100: if (flg) {
101: PetscViewerASCIIOpen(eps->comm,monfilename,&monviewer);
102: EPSSetMonitor(eps,EPSDefaultMonitor,monviewer);
103: }
104: PetscOptionsName("-eps_xmonitor","Monitor error estimates graphically","EPSSetMonitor",&flg);
105: if (flg) {
106: EPSSetMonitor(eps,EPSLGMonitor,PETSC_NULL);
107: }
108: /* -----------------------------------------------------------------------*/
109: PetscOptionsTruthGroupBegin("-eps_largest_magnitude","compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
110: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE);}
111: PetscOptionsTruthGroup("-eps_smallest_magnitude","compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
112: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE);}
113: PetscOptionsTruthGroup("-eps_largest_real","compute largest real parts","EPSSetWhichEigenpairs",&flg);
114: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL);}
115: PetscOptionsTruthGroup("-eps_smallest_real","compute smallest real parts","EPSSetWhichEigenpairs",&flg);
116: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);}
117: PetscOptionsTruthGroup("-eps_largest_imaginary","compute largest imaginary parts","EPSSetWhichEigenpairs",&flg);
118: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY);}
119: PetscOptionsTruthGroupEnd("-eps_smallest_imaginary","compute smallest imaginary parts","EPSSetWhichEigenpairs",&flg);
120: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY);}
122: PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",0);
123: PetscOptionsName("-eps_view_binary","Save the matrices associated to the eigenproblem","EPSSetFromOptions",0);
124: PetscOptionsName("-eps_plot_eigs","Make a plot of the computed eigenvalues","EPSSolve",0);
125:
126: if (eps->ops->setfromoptions) {
127: (*eps->ops->setfromoptions)(eps);
128: }
129: PetscOptionsEnd();
131: STSetFromOptions(eps->OP);
132: return(0);
133: }
137: /*@
138: EPSGetTolerances - Gets the tolerance and maximum
139: iteration count used by the default EPS convergence tests.
141: Not Collective
143: Input Parameter:
144: . eps - the eigensolver context
145:
146: Output Parameters:
147: + tol - the convergence tolerance
148: - maxits - maximum number of iterations
150: Notes:
151: The user can specify PETSC_NULL for any parameter that is not needed.
153: Level: intermediate
155: .seealso: EPSSetTolerances()
156: @*/
157: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,int *maxits)
158: {
161: if (tol) *tol = eps->tol;
162: if (maxits) *maxits = eps->max_it;
163: return(0);
164: }
168: /*@
169: EPSSetTolerances - Sets the tolerance and maximum
170: iteration count used by the default EPS convergence testers.
172: Collective on EPS
174: Input Parameters:
175: + eps - the eigensolver context
176: . tol - the convergence tolerance
177: - maxits - maximum number of iterations to use
179: Options Database Keys:
180: + -eps_tol <tol> - Sets the convergence tolerance
181: - -eps_max_it <maxits> - Sets the maximum number of iterations allowed
183: Notes:
184: Use PETSC_DEFAULT to retain the default value of any of the tolerances.
186: Level: intermediate
188: .seealso: EPSGetTolerances()
189: @*/
190: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,int maxits)
191: {
194: if (tol != PETSC_DEFAULT) eps->tol = tol;
195: if (maxits != PETSC_DEFAULT) eps->max_it = maxits;
196: return(0);
197: }
201: /*@
202: EPSGetDimensions - Gets the number of eigenvalues to compute
203: and the dimension of the subspace.
205: Not Collective
207: Input Parameter:
208: . eps - the eigensolver context
209:
210: Output Parameters:
211: + nev - number of eigenvalues to compute
212: - ncv - the maximum dimension of the subspace to be used by the solver
214: Notes:
215: The user can specify PETSC_NULL for any parameter that is not needed.
217: Level: intermediate
219: .seealso: EPSSetDimensions()
220: @*/
221: PetscErrorCode EPSGetDimensions(EPS eps,int *nev,int *ncv)
222: {
225: if( nev ) *nev = eps->nev;
226: if( ncv ) *ncv = eps->ncv;
227: return(0);
228: }
232: /*@
233: EPSSetDimensions - Sets the number of eigenvalues to compute
234: and the dimension of the subspace.
236: Collective on EPS
238: Input Parameters:
239: + eps - the eigensolver context
240: . nev - number of eigenvalues to compute
241: - ncv - the maximum dimension of the subspace to be used by the solver
243: Options Database Keys:
244: + -eps_nev <nev> - Sets the number of eigenvalues
245: - -eps_ncv <ncv> - Sets the dimension of the subspace
247: Notes:
248: Use PETSC_DEFAULT to retain the previous value of any parameter.
250: Use PETSC_DECIDE for ncv to assign a reasonably good value, which is
251: dependent on the solution method.
253: Level: intermediate
255: .seealso: EPSGetDimensions()
256: @*/
257: PetscErrorCode EPSSetDimensions(EPS eps,int nev,int ncv)
258: {
262: if( nev != PETSC_DEFAULT ) {
263: if (nev<1) SETERRQ(1,"Illegal value of nev. Must be > 0");
264: eps->nev = nev;
265: eps->setupcalled = 0;
266: }
267: if( ncv != PETSC_DEFAULT ) {
268: if( ncv == PETSC_DECIDE ) eps->ncv = 0;
269: else {
270: if (ncv<1) SETERRQ(1,"Illegal value of ncv. Must be > 0");
271: eps->ncv = ncv;
272: }
273: eps->setupcalled = 0;
274: }
275: return(0);
276: }
280: /*@
281: EPSSetWhichEigenpairs - Specifies which portion of the spectrum is
282: to be sought.
284: Collective on EPS
286: Input Parameter:
287: . eps - eigensolver context obtained from EPSCreate()
289: Output Parameter:
290: . which - the portion of the spectrum to be sought
292: Possible values:
293: The parameter 'which' can have one of these values:
294:
295: + EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
296: . EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
297: . EPS_LARGEST_REAL - largest real parts
298: . EPS_SMALLEST_REAL - smallest real parts
299: . EPS_LARGEST_IMAGINARY - largest imaginary parts
300: - EPS_SMALLEST_IMAGINARY - smallest imaginary parts
302: Options Database Keys:
303: + -eps_largest_magnitude - Sets largest eigenvalues in magnitude
304: . -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
305: . -eps_largest_real - Sets largest real parts
306: . -eps_smallest_real - Sets smallest real parts
307: . -eps_largest_imaginary - Sets largest imaginary parts in magnitude
308: - -eps_smallest_imaginary - Sets smallest imaginary parts in magnitude
310: Notes:
311: Not all eigensolvers implemented in EPS account for all the possible values
312: stated above. Also, some values make sense only for certain types of
313: problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY
314: and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part
315: for eigenvalue selection.
316:
317: Level: intermediate
319: .seealso: EPSGetWhichEigenpairs(), EPSSortEigenvalues()
320: @*/
321: PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)
322: {
325: eps->which = which;
326: return(0);
327: }
331: /*@C
332: EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be
333: sought.
335: Not Collective
337: Input Parameter:
338: . eps - eigensolver context obtained from EPSCreate()
340: Output Parameter:
341: . which - the portion of the spectrum to be sought
343: Notes:
344: See EPSSetWhichEigenpairs() for possible values of which
346: Level: intermediate
348: .seealso: EPSSetWhichEigenpairs()
349: @*/
350: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)
351: {
354: *which = eps->which;
355: return(0);
356: }
360: /*@
361: EPSSetProblemType - Specifies the type of the eigenvalue problem.
363: Collective on EPS
365: Input Parameters:
366: + eps - the eigensolver context
367: - type - a known type of eigenvalue problem
369: Options Database Keys:
370: + -eps_hermitian - Hermitian eigenvalue problem
371: . -eps_gen_hermitian - generalized Hermitian eigenvalue problem
372: . -eps_non_hermitian - non-Hermitian eigenvalue problem
373: - -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
374:
375: Note:
376: Allowed values for the problem type are: Hermitian (EPS_HEP), non-Hermitian
377: (EPS_NHEP), generalized Hermitian (EPS_GHEP) and generalized non-Hermitian
378: (EPS_GNHEP).
380: This function must be used to instruct SLEPc to exploit symmetry. If no
381: problem type is specified, by default a non-Hermitian problem is assumed
382: (either standard or generalized). If the user knows that the problem is
383: Hermitian (i.e. A=A^H) of generalized Hermitian (i.e. A=A^H, B=B^H, and
384: B positive definite) then it is recommended to set the problem type so
385: that eigensolver can exploit these properties.
387: Level: beginner
389: .seealso: EPSSetOperators(), EPSSetType(), EPSProblemType
390: @*/
391: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)
392: {
394: Mat A,B;
399: STGetOperators(eps->OP,&A,&B);
400: if (!A) { SETERRQ(1,"Must call EPSSetOperators() first"); }
402: switch (type) {
403: case EPS_HEP:
404: eps->isgeneralized = PETSC_FALSE;
405: eps->ishermitian = PETSC_TRUE;
406: STSetBilinearForm(eps->OP,STINNER_HERMITIAN);
407: break;
408: case EPS_NHEP:
409: eps->isgeneralized = PETSC_FALSE;
410: eps->ishermitian = PETSC_FALSE;
411: STSetBilinearForm(eps->OP,STINNER_HERMITIAN);
412: break;
413: case EPS_GHEP:
414: eps->isgeneralized = PETSC_TRUE;
415: eps->ishermitian = PETSC_TRUE;
416: STSetBilinearForm(eps->OP,STINNER_B_HERMITIAN);
417: break;
418: case EPS_GNHEP:
419: eps->isgeneralized = PETSC_TRUE;
420: eps->ishermitian = PETSC_FALSE;
421: STSetBilinearForm(eps->OP,STINNER_HERMITIAN);
422: break;
423: /*
424: case EPS_CSEP:
425: eps->isgeneralized = PETSC_FALSE;
426: eps->ishermitian = PETSC_FALSE;
427: STSetBilinearForm(eps->OP,STINNER_SYMMETRIC);
428: break;
429: case EPS_GCSEP:
430: eps->isgeneralized = PETSC_TRUE;
431: eps->ishermitian = PETSC_FALSE;
432: STSetBilinearForm(eps->OP,STINNER_B_SYMMETRIC);
433: break;
434: */
435: default:
436: SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
437: }
438: eps->problem_type = type;
440: return(0);
441: }
445: /*@C
446: EPSGetProblemType - Gets the problem type from the EPS object.
448: Not Collective
450: Input Parameter:
451: . eps - the eigensolver context
453: Output Parameter:
454: . type - name of EPS problem type
456: Level: intermediate
458: .seealso: EPSSetProblemType()
459: @*/
460: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)
461: {
464: *type = eps->problem_type;
465: return(0);
466: }
470: /*@
471: EPSSetClass - Specifies the eigensolver class: either one-sided or two-sided.
473: Collective on EPS
475: Input Parameters:
476: + eps - the eigensolver context
477: - class - the class of solver
479: Options Database Keys:
480: + -eps_oneside - one-sided solver
481: - -eps_twoside - two-sided solver
482:
483: Note:
484: Allowed solver classes are: one-sided (EPS_ONE_SIDE) and two-sided (EPS_TWO_SIDE).
485: One-sided eigensolvers are the standard ones, which allow the computation of
486: eigenvalues and (right) eigenvectors, whereas two-sided eigensolvers compute
487: left eigenvectors as well.
489: Level: beginner
491: .seealso: EPSGetLeftVector(), EPSComputeRelativeErrorLeft(), EPSSetLeftInitialVector(),
492: EPSClass
493: @*/
494: PetscErrorCode EPSSetClass(EPS eps,EPSClass cl)
495: {
501: if (cl != EPS_ONE_SIDE && cl != EPS_TWO_SIDE) SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown eigensolver class");
502: if (eps->solverclass!=cl) {
503: if (eps->solverclass == EPS_TWO_SIDE) { EPSFreeSolution(eps); }
504: eps->solverclass = cl;
505: }
507: return(0);
508: }
512: /*@C
513: EPSGetClass - Gets the eigensolver class from the EPS object.
515: Not Collective
517: Input Parameter:
518: . eps - the eigensolver context
520: Output Parameter:
521: . class - class of EPS solver (either one-sided or two-sided)
523: Level: intermediate
525: .seealso: EPSSetClass()
526: @*/
527: PetscErrorCode EPSGetClass(EPS eps,EPSClass *cl)
528: {
531: *cl = eps->solverclass;
532: return(0);
533: }
537: /*@
538: EPSSetOrthogonalization - Specifies the type of orthogonalization technique
539: to be used inside the eigensolver (classical or modified Gram-Schmidt with
540: or without refinement).
542: Collective on EPS
544: Input Parameters:
545: + eps - the eigensolver context
546: . type - the type of orthogonalization technique
547: . refinement - type of refinement
548: - eta - parameter for selective refinement
550: Options Database Keys:
551: + -eps_orthog_type <type> - Where <type> is cgs for Classical Gram-Schmidt
552: orthogonalization (default)
553: or mgs for Modified Gram-Schmidt orthogonalization
554: . -eps_orthog_refinement <type> - Where <type> is one of never, ifneeded
555: (default) or always
556: - -eps_orthog_eta <eta> - For setting the value of eta
557:
558: Notes:
559: The default settings work well for most problems.
561: The parameter eta should be a real value between 0 and 1 (or PETSC_DEFAULT).
562: The value of eta is used only when the refinement type is "ifneeded".
564: When using several processors, MGS is likely to result in bad scalability.
566: Level: advanced
568: .seealso: EPSOrthogonalize(), EPSGetOrthogonalization()
569: @*/
570: PetscErrorCode EPSSetOrthogonalization(EPS eps,EPSOrthogonalizationType type, EPSOrthogonalizationRefinementType refinement, PetscReal eta)
571: {
574: switch (type) {
575: case EPS_CGS_ORTH:
576: case EPS_MGS_ORTH:
577: eps->orthog_type = type;
578: break;
579: default:
580: SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown orthogonalization type");
581: }
582: switch (refinement) {
583: case EPS_ORTH_REFINE_NEVER:
584: case EPS_ORTH_REFINE_IFNEEDED:
585: case EPS_ORTH_REFINE_ALWAYS:
586: eps->orthog_ref = refinement;
587: break;
588: default:
589: SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown refinement type");
590: }
591: if (eta != PETSC_DEFAULT && eta <= 0.0 && eta > 1.0) {
592: SETERRQ(PETSC_ERR_ARG_WRONG,"Invalid eta value");
593: }
594: eps->orthog_eta = eta;
595: return(0);
596: }
600: /*@C
601: EPSGetOrthogonalization - Gets the orthogonalization settings from the
602: EPS object.
604: Not Collective
606: Input Parameter:
607: . eps - Eigensolver context
609: Output Parameter:
610: + type - type of orthogonalization technique
611: . refinement - type of refinement
612: - eta - parameter for selective refinement
614: Level: advanced
616: .seealso: EPSOrthogonalize(), EPSSetOrthogonalization()
617: @*/
618: PetscErrorCode EPSGetOrthogonalization(EPS eps,EPSOrthogonalizationType *type,EPSOrthogonalizationRefinementType *refinement, PetscReal *eta)
619: {
622: if (type) *type = eps->orthog_type;
623: if (refinement) *refinement = eps->orthog_ref;
624: if (eta) *eta = eps->orthog_eta;
625: return(0);
626: }
630: /*@C
631: EPSSetOptionsPrefix - Sets the prefix used for searching for all
632: EPS options in the database.
634: Collective on EPS
636: Input Parameters:
637: + eps - the eigensolver context
638: - prefix - the prefix string to prepend to all EPS option requests
640: Notes:
641: A hyphen (-) must NOT be given at the beginning of the prefix name.
642: The first character of all runtime options is AUTOMATICALLY the
643: hyphen.
645: For example, to distinguish between the runtime options for two
646: different EPS contexts, one could call
647: .vb
648: EPSSetOptionsPrefix(eps1,"eig1_")
649: EPSSetOptionsPrefix(eps2,"eig2_")
650: .ve
652: Level: advanced
654: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
655: @*/
656: PetscErrorCode EPSSetOptionsPrefix(EPS eps,char *prefix)
657: {
661: PetscObjectSetOptionsPrefix((PetscObject)eps, prefix);
662: STSetOptionsPrefix(eps->OP,prefix);
663: return(0);
664: }
665:
668: /*@C
669: EPSAppendOptionsPrefix - Appends to the prefix used for searching for all
670: EPS options in the database.
672: Collective on EPS
674: Input Parameters:
675: + eps - the eigensolver context
676: - prefix - the prefix string to prepend to all EPS option requests
678: Notes:
679: A hyphen (-) must NOT be given at the beginning of the prefix name.
680: The first character of all runtime options is AUTOMATICALLY the hyphen.
682: Level: advanced
684: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
685: @*/
686: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,char *prefix)
687: {
691: PetscObjectAppendOptionsPrefix((PetscObject)eps, prefix);
692: STAppendOptionsPrefix(eps->OP,prefix);
693: return(0);
694: }
698: /*@C
699: EPSGetOptionsPrefix - Gets the prefix used for searching for all
700: EPS options in the database.
702: Not Collective
704: Input Parameters:
705: . eps - the eigensolver context
707: Output Parameters:
708: . prefix - pointer to the prefix string used is returned
710: Notes: On the fortran side, the user should pass in a string 'prefix' of
711: sufficient length to hold the prefix.
713: Level: advanced
715: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
716: @*/
717: PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])
718: {
722: PetscObjectGetOptionsPrefix((PetscObject)eps, prefix);
723: return(0);
724: }