Line data Source code
1 : /*
2 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3 : SLEPc - Scalable Library for Eigenvalue Problem Computations
4 : Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5 :
6 : This file is part of SLEPc.
7 : SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9 : */
10 : /*
11 : EPS routines related to options that can be set via the command-line
12 : or procedurally.
13 : */
14 :
15 : #include <slepc/private/epsimpl.h> /*I "slepceps.h" I*/
16 : #include <petscdraw.h>
17 :
18 : /*@C
19 : EPSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
20 : indicated by the user.
21 :
22 : Collective
23 :
24 : Input Parameters:
25 : + eps - the eigensolver context
26 : . opt - the command line option for this monitor
27 : . name - the monitor type one is seeking
28 : . ctx - an optional user context for the monitor, or NULL
29 : - trackall - whether this monitor tracks all eigenvalues or not
30 :
31 : Level: developer
32 :
33 : .seealso: EPSMonitorSet(), EPSSetTrackAll()
34 : @*/
35 1641 : PetscErrorCode EPSMonitorSetFromOptions(EPS eps,const char opt[],const char name[],void *ctx,PetscBool trackall)
36 : {
37 1641 : PetscErrorCode (*mfunc)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
38 1641 : PetscErrorCode (*cfunc)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**);
39 1641 : PetscErrorCode (*dfunc)(PetscViewerAndFormat**);
40 1641 : PetscViewerAndFormat *vf;
41 1641 : PetscViewer viewer;
42 1641 : PetscViewerFormat format;
43 1641 : PetscViewerType vtype;
44 1641 : char key[PETSC_MAX_PATH_LEN];
45 1641 : PetscBool flg;
46 :
47 1641 : PetscFunctionBegin;
48 1641 : PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,opt,&viewer,&format,&flg));
49 1641 : if (!flg) PetscFunctionReturn(PETSC_SUCCESS);
50 :
51 8 : PetscCall(PetscViewerGetType(viewer,&vtype));
52 8 : PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
53 8 : PetscCall(PetscFunctionListFind(EPSMonitorList,key,&mfunc));
54 8 : PetscCheck(mfunc,PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Specified viewer and format not supported");
55 8 : PetscCall(PetscFunctionListFind(EPSMonitorCreateList,key,&cfunc));
56 8 : PetscCall(PetscFunctionListFind(EPSMonitorDestroyList,key,&dfunc));
57 8 : if (!cfunc) cfunc = PetscViewerAndFormatCreate_Internal;
58 8 : if (!dfunc) dfunc = PetscViewerAndFormatDestroy;
59 :
60 8 : PetscCall((*cfunc)(viewer,format,ctx,&vf));
61 8 : PetscCall(PetscOptionsRestoreViewer(&viewer));
62 8 : PetscCall(EPSMonitorSet(eps,mfunc,vf,(PetscErrorCode(*)(void **))dfunc));
63 8 : if (trackall) PetscCall(EPSSetTrackAll(eps,PETSC_TRUE));
64 8 : PetscFunctionReturn(PETSC_SUCCESS);
65 : }
66 :
67 : /*@
68 : EPSSetFromOptions - Sets EPS options from the options database.
69 : This routine must be called before EPSSetUp() if the user is to be
70 : allowed to set the solver type.
71 :
72 : Collective
73 :
74 : Input Parameters:
75 : . eps - the eigensolver context
76 :
77 : Notes:
78 : To see all options, run your program with the -help option.
79 :
80 : Level: beginner
81 :
82 : .seealso: EPSSetOptionsPrefix()
83 : @*/
84 547 : PetscErrorCode EPSSetFromOptions(EPS eps)
85 : {
86 547 : char type[256];
87 547 : PetscBool set,flg,flg1,flg2,flg3,bval;
88 547 : PetscReal r,array[2]={0,0};
89 547 : PetscScalar s;
90 547 : PetscInt i,j,k;
91 547 : EPSBalance bal;
92 :
93 547 : PetscFunctionBegin;
94 547 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
95 547 : PetscCall(EPSRegisterAll());
96 1641 : PetscObjectOptionsBegin((PetscObject)eps);
97 603 : PetscCall(PetscOptionsFList("-eps_type","Eigensolver method","EPSSetType",EPSList,(char*)(((PetscObject)eps)->type_name?((PetscObject)eps)->type_name:EPSKRYLOVSCHUR),type,sizeof(type),&flg));
98 547 : if (flg) PetscCall(EPSSetType(eps,type));
99 295 : else if (!((PetscObject)eps)->type_name) PetscCall(EPSSetType(eps,EPSKRYLOVSCHUR));
100 :
101 547 : PetscCall(PetscOptionsBoolGroupBegin("-eps_hermitian","Hermitian eigenvalue problem","EPSSetProblemType",&flg));
102 547 : if (flg) PetscCall(EPSSetProblemType(eps,EPS_HEP));
103 547 : PetscCall(PetscOptionsBoolGroup("-eps_gen_hermitian","Generalized Hermitian eigenvalue problem","EPSSetProblemType",&flg));
104 547 : if (flg) PetscCall(EPSSetProblemType(eps,EPS_GHEP));
105 547 : PetscCall(PetscOptionsBoolGroup("-eps_non_hermitian","Non-Hermitian eigenvalue problem","EPSSetProblemType",&flg));
106 547 : if (flg) PetscCall(EPSSetProblemType(eps,EPS_NHEP));
107 547 : PetscCall(PetscOptionsBoolGroup("-eps_gen_non_hermitian","Generalized non-Hermitian eigenvalue problem","EPSSetProblemType",&flg));
108 547 : if (flg) PetscCall(EPSSetProblemType(eps,EPS_GNHEP));
109 547 : PetscCall(PetscOptionsBoolGroup("-eps_pos_gen_non_hermitian","Generalized non-Hermitian eigenvalue problem with positive semi-definite B","EPSSetProblemType",&flg));
110 547 : if (flg) PetscCall(EPSSetProblemType(eps,EPS_PGNHEP));
111 547 : PetscCall(PetscOptionsBoolGroupEnd("-eps_gen_indefinite","Generalized Hermitian-indefinite eigenvalue problem","EPSSetProblemType",&flg));
112 547 : if (flg) PetscCall(EPSSetProblemType(eps,EPS_GHIEP));
113 :
114 547 : PetscCall(PetscOptionsBoolGroupBegin("-eps_ritz","Rayleigh-Ritz extraction","EPSSetExtraction",&flg));
115 547 : if (flg) PetscCall(EPSSetExtraction(eps,EPS_RITZ));
116 547 : PetscCall(PetscOptionsBoolGroup("-eps_harmonic","Harmonic Ritz extraction","EPSSetExtraction",&flg));
117 547 : if (flg) PetscCall(EPSSetExtraction(eps,EPS_HARMONIC));
118 547 : PetscCall(PetscOptionsBoolGroup("-eps_harmonic_relative","Relative harmonic Ritz extraction","EPSSetExtraction",&flg));
119 547 : if (flg) PetscCall(EPSSetExtraction(eps,EPS_HARMONIC_RELATIVE));
120 547 : PetscCall(PetscOptionsBoolGroup("-eps_harmonic_right","Right harmonic Ritz extraction","EPSSetExtraction",&flg));
121 547 : if (flg) PetscCall(EPSSetExtraction(eps,EPS_HARMONIC_RIGHT));
122 547 : PetscCall(PetscOptionsBoolGroup("-eps_harmonic_largest","Largest harmonic Ritz extraction","EPSSetExtraction",&flg));
123 547 : if (flg) PetscCall(EPSSetExtraction(eps,EPS_HARMONIC_LARGEST));
124 547 : PetscCall(PetscOptionsBoolGroup("-eps_refined","Refined Ritz extraction","EPSSetExtraction",&flg));
125 547 : if (flg) PetscCall(EPSSetExtraction(eps,EPS_REFINED));
126 547 : PetscCall(PetscOptionsBoolGroupEnd("-eps_refined_harmonic","Refined harmonic Ritz extraction","EPSSetExtraction",&flg));
127 547 : if (flg) PetscCall(EPSSetExtraction(eps,EPS_REFINED_HARMONIC));
128 :
129 547 : bal = eps->balance;
130 547 : PetscCall(PetscOptionsEnum("-eps_balance","Balancing method","EPSSetBalance",EPSBalanceTypes,(PetscEnum)bal,(PetscEnum*)&bal,&flg1));
131 547 : j = eps->balance_its;
132 547 : PetscCall(PetscOptionsInt("-eps_balance_its","Number of iterations in balancing","EPSSetBalance",eps->balance_its,&j,&flg2));
133 547 : r = eps->balance_cutoff;
134 547 : PetscCall(PetscOptionsReal("-eps_balance_cutoff","Cutoff value in balancing","EPSSetBalance",eps->balance_cutoff,&r,&flg3));
135 547 : if (flg1 || flg2 || flg3) PetscCall(EPSSetBalance(eps,bal,j,r));
136 :
137 547 : i = eps->max_it;
138 547 : PetscCall(PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&i,&flg1));
139 547 : r = eps->tol;
140 897 : PetscCall(PetscOptionsReal("-eps_tol","Tolerance","EPSSetTolerances",SlepcDefaultTol(eps->tol),&r,&flg2));
141 547 : if (flg1 || flg2) PetscCall(EPSSetTolerances(eps,r,i));
142 :
143 547 : PetscCall(PetscOptionsBoolGroupBegin("-eps_conv_rel","Relative error convergence test","EPSSetConvergenceTest",&flg));
144 547 : if (flg) PetscCall(EPSSetConvergenceTest(eps,EPS_CONV_REL));
145 547 : PetscCall(PetscOptionsBoolGroup("-eps_conv_norm","Convergence test relative to the eigenvalue and the matrix norms","EPSSetConvergenceTest",&flg));
146 547 : if (flg) PetscCall(EPSSetConvergenceTest(eps,EPS_CONV_NORM));
147 547 : PetscCall(PetscOptionsBoolGroup("-eps_conv_abs","Absolute error convergence test","EPSSetConvergenceTest",&flg));
148 547 : if (flg) PetscCall(EPSSetConvergenceTest(eps,EPS_CONV_ABS));
149 547 : PetscCall(PetscOptionsBoolGroupEnd("-eps_conv_user","User-defined convergence test","EPSSetConvergenceTest",&flg));
150 547 : if (flg) PetscCall(EPSSetConvergenceTest(eps,EPS_CONV_USER));
151 :
152 547 : PetscCall(PetscOptionsBoolGroupBegin("-eps_stop_basic","Stop iteration if all eigenvalues converged or max_it reached","EPSSetStoppingTest",&flg));
153 547 : if (flg) PetscCall(EPSSetStoppingTest(eps,EPS_STOP_BASIC));
154 547 : PetscCall(PetscOptionsBoolGroupEnd("-eps_stop_user","User-defined stopping test","EPSSetStoppingTest",&flg));
155 547 : if (flg) PetscCall(EPSSetStoppingTest(eps,EPS_STOP_USER));
156 :
157 547 : i = eps->nev;
158 547 : PetscCall(PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&i,&flg1));
159 547 : j = eps->ncv;
160 547 : PetscCall(PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&j,&flg2));
161 547 : k = eps->mpd;
162 547 : PetscCall(PetscOptionsInt("-eps_mpd","Maximum dimension of projected problem","EPSSetDimensions",eps->mpd,&k,&flg3));
163 547 : if (flg1 || flg2 || flg3) PetscCall(EPSSetDimensions(eps,i,j,k));
164 :
165 547 : PetscCall(PetscOptionsBoolGroupBegin("-eps_largest_magnitude","Compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg));
166 547 : if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE));
167 547 : PetscCall(PetscOptionsBoolGroup("-eps_smallest_magnitude","Compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg));
168 547 : if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE));
169 547 : PetscCall(PetscOptionsBoolGroup("-eps_largest_real","Compute eigenvalues with largest real parts","EPSSetWhichEigenpairs",&flg));
170 547 : if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL));
171 547 : PetscCall(PetscOptionsBoolGroup("-eps_smallest_real","Compute eigenvalues with smallest real parts","EPSSetWhichEigenpairs",&flg));
172 547 : if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL));
173 547 : PetscCall(PetscOptionsBoolGroup("-eps_largest_imaginary","Compute eigenvalues with largest imaginary parts","EPSSetWhichEigenpairs",&flg));
174 547 : if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY));
175 547 : PetscCall(PetscOptionsBoolGroup("-eps_smallest_imaginary","Compute eigenvalues with smallest imaginary parts","EPSSetWhichEigenpairs",&flg));
176 547 : if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY));
177 547 : PetscCall(PetscOptionsBoolGroup("-eps_target_magnitude","Compute eigenvalues closest to target","EPSSetWhichEigenpairs",&flg));
178 547 : if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE));
179 547 : PetscCall(PetscOptionsBoolGroup("-eps_target_real","Compute eigenvalues with real parts closest to target","EPSSetWhichEigenpairs",&flg));
180 547 : if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_REAL));
181 547 : PetscCall(PetscOptionsBoolGroup("-eps_target_imaginary","Compute eigenvalues with imaginary parts closest to target","EPSSetWhichEigenpairs",&flg));
182 547 : if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_IMAGINARY));
183 547 : PetscCall(PetscOptionsBoolGroupEnd("-eps_all","Compute all eigenvalues in an interval or a region","EPSSetWhichEigenpairs",&flg));
184 547 : if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_ALL));
185 :
186 547 : PetscCall(PetscOptionsScalar("-eps_target","Value of the target","EPSSetTarget",eps->target,&s,&flg));
187 547 : if (flg) {
188 35 : if (eps->which!=EPS_TARGET_REAL && eps->which!=EPS_TARGET_IMAGINARY) PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE));
189 35 : PetscCall(EPSSetTarget(eps,s));
190 : }
191 :
192 547 : k = 2;
193 547 : PetscCall(PetscOptionsRealArray("-eps_interval","Computational interval (two real values separated with a comma without spaces)","EPSSetInterval",array,&k,&flg));
194 547 : if (flg) {
195 11 : PetscCheck(k>1,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_SIZ,"Must pass two values in -eps_interval (comma-separated without spaces)");
196 11 : PetscCall(EPSSetWhichEigenpairs(eps,EPS_ALL));
197 11 : PetscCall(EPSSetInterval(eps,array[0],array[1]));
198 : }
199 :
200 547 : PetscCall(PetscOptionsBool("-eps_true_residual","Compute true residuals explicitly","EPSSetTrueResidual",eps->trueres,&eps->trueres,NULL));
201 547 : PetscCall(PetscOptionsBool("-eps_purify","Postprocess eigenvectors for purification","EPSSetPurify",eps->purify,&bval,&flg));
202 547 : if (flg) PetscCall(EPSSetPurify(eps,bval));
203 547 : PetscCall(PetscOptionsBool("-eps_two_sided","Use two-sided variant (to compute left eigenvectors)","EPSSetTwoSided",eps->twosided,&bval,&flg));
204 547 : if (flg) PetscCall(EPSSetTwoSided(eps,bval));
205 :
206 : /* -----------------------------------------------------------------------*/
207 : /*
208 : Cancels all monitors hardwired into code before call to EPSSetFromOptions()
209 : */
210 547 : PetscCall(PetscOptionsBool("-eps_monitor_cancel","Remove any hardwired monitor routines","EPSMonitorCancel",PETSC_FALSE,&flg,&set));
211 547 : if (set && flg) PetscCall(EPSMonitorCancel(eps));
212 547 : PetscCall(EPSMonitorSetFromOptions(eps,"-eps_monitor","first_approximation",NULL,PETSC_FALSE));
213 547 : PetscCall(EPSMonitorSetFromOptions(eps,"-eps_monitor_all","all_approximations",NULL,PETSC_TRUE));
214 547 : PetscCall(EPSMonitorSetFromOptions(eps,"-eps_monitor_conv","convergence_history",NULL,PETSC_FALSE));
215 :
216 : /* -----------------------------------------------------------------------*/
217 547 : PetscCall(PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",&set));
218 547 : PetscCall(PetscOptionsName("-eps_view_vectors","View computed eigenvectors","EPSVectorsView",&set));
219 547 : PetscCall(PetscOptionsName("-eps_view_values","View computed eigenvalues","EPSValuesView",&set));
220 547 : PetscCall(PetscOptionsName("-eps_converged_reason","Print reason for convergence, and number of iterations","EPSConvergedReasonView",&set));
221 547 : PetscCall(PetscOptionsName("-eps_error_absolute","Print absolute errors of each eigenpair","EPSErrorView",&set));
222 547 : PetscCall(PetscOptionsName("-eps_error_relative","Print relative errors of each eigenpair","EPSErrorView",&set));
223 547 : PetscCall(PetscOptionsName("-eps_error_backward","Print backward errors of each eigenpair","EPSErrorView",&set));
224 :
225 547 : PetscTryTypeMethod(eps,setfromoptions,PetscOptionsObject);
226 547 : PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)eps,PetscOptionsObject));
227 547 : PetscOptionsEnd();
228 :
229 547 : if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
230 547 : PetscCall(BVSetFromOptions(eps->V));
231 547 : if (!eps->rg) PetscCall(EPSGetRG(eps,&eps->rg));
232 547 : PetscCall(RGSetFromOptions(eps->rg));
233 547 : if (eps->useds) {
234 499 : if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
235 499 : PetscCall(EPSSetDSType(eps));
236 499 : PetscCall(DSSetFromOptions(eps->ds));
237 : }
238 547 : if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
239 547 : PetscCall(EPSSetDefaultST(eps));
240 547 : PetscCall(STSetFromOptions(eps->st));
241 547 : PetscFunctionReturn(PETSC_SUCCESS);
242 : }
243 :
244 : /*@C
245 : EPSGetTolerances - Gets the tolerance and maximum iteration count used
246 : by the EPS convergence tests.
247 :
248 : Not Collective
249 :
250 : Input Parameter:
251 : . eps - the eigensolver context
252 :
253 : Output Parameters:
254 : + tol - the convergence tolerance
255 : - maxits - maximum number of iterations
256 :
257 : Notes:
258 : The user can specify NULL for any parameter that is not needed.
259 :
260 : Level: intermediate
261 :
262 : .seealso: EPSSetTolerances()
263 : @*/
264 251 : PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,PetscInt *maxits)
265 : {
266 251 : PetscFunctionBegin;
267 251 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
268 251 : if (tol) *tol = eps->tol;
269 251 : if (maxits) *maxits = eps->max_it;
270 251 : PetscFunctionReturn(PETSC_SUCCESS);
271 : }
272 :
273 : /*@
274 : EPSSetTolerances - Sets the tolerance and maximum iteration count used
275 : by the EPS convergence tests.
276 :
277 : Logically Collective
278 :
279 : Input Parameters:
280 : + eps - the eigensolver context
281 : . tol - the convergence tolerance
282 : - maxits - maximum number of iterations to use
283 :
284 : Options Database Keys:
285 : + -eps_tol <tol> - Sets the convergence tolerance
286 : - -eps_max_it <maxits> - Sets the maximum number of iterations allowed
287 :
288 : Notes:
289 : Use PETSC_DEFAULT for either argument to assign a reasonably good value.
290 :
291 : Level: intermediate
292 :
293 : .seealso: EPSGetTolerances()
294 : @*/
295 525 : PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,PetscInt maxits)
296 : {
297 525 : PetscFunctionBegin;
298 525 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
299 2100 : PetscValidLogicalCollectiveReal(eps,tol,2);
300 2100 : PetscValidLogicalCollectiveInt(eps,maxits,3);
301 525 : if (tol == (PetscReal)PETSC_DEFAULT) {
302 19 : eps->tol = PETSC_DEFAULT;
303 19 : eps->state = EPS_STATE_INITIAL;
304 : } else {
305 506 : PetscCheck(tol>0.0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
306 506 : eps->tol = tol;
307 : }
308 525 : if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
309 322 : eps->max_it = PETSC_DEFAULT;
310 322 : eps->state = EPS_STATE_INITIAL;
311 : } else {
312 203 : PetscCheck(maxits>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
313 203 : eps->max_it = maxits;
314 : }
315 525 : PetscFunctionReturn(PETSC_SUCCESS);
316 : }
317 :
318 : /*@C
319 : EPSGetDimensions - Gets the number of eigenvalues to compute
320 : and the dimension of the subspace.
321 :
322 : Not Collective
323 :
324 : Input Parameter:
325 : . eps - the eigensolver context
326 :
327 : Output Parameters:
328 : + nev - number of eigenvalues to compute
329 : . ncv - the maximum dimension of the subspace to be used by the solver
330 : - mpd - the maximum dimension allowed for the projected problem
331 :
332 : Level: intermediate
333 :
334 : .seealso: EPSSetDimensions()
335 : @*/
336 368 : PetscErrorCode EPSGetDimensions(EPS eps,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
337 : {
338 368 : PetscFunctionBegin;
339 368 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
340 368 : if (nev) *nev = eps->nev;
341 368 : if (ncv) *ncv = eps->ncv;
342 368 : if (mpd) *mpd = eps->mpd;
343 368 : PetscFunctionReturn(PETSC_SUCCESS);
344 : }
345 :
346 : /*@
347 : EPSSetDimensions - Sets the number of eigenvalues to compute
348 : and the dimension of the subspace.
349 :
350 : Logically Collective
351 :
352 : Input Parameters:
353 : + eps - the eigensolver context
354 : . nev - number of eigenvalues to compute
355 : . ncv - the maximum dimension of the subspace to be used by the solver
356 : - mpd - the maximum dimension allowed for the projected problem
357 :
358 : Options Database Keys:
359 : + -eps_nev <nev> - Sets the number of eigenvalues
360 : . -eps_ncv <ncv> - Sets the dimension of the subspace
361 : - -eps_mpd <mpd> - Sets the maximum projected dimension
362 :
363 : Notes:
364 : Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
365 : dependent on the solution method.
366 :
367 : The parameters ncv and mpd are intimately related, so that the user is advised
368 : to set one of them at most. Normal usage is that
369 : (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
370 : (b) in cases where nev is large, the user sets mpd.
371 :
372 : The value of ncv should always be between nev and (nev+mpd), typically
373 : ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
374 : a smaller value should be used.
375 :
376 : When computing all eigenvalues in an interval, see EPSSetInterval(), these
377 : parameters lose relevance, and tuning must be done with
378 : EPSKrylovSchurSetDimensions().
379 :
380 : Level: intermediate
381 :
382 : .seealso: EPSGetDimensions(), EPSSetInterval(), EPSKrylovSchurSetDimensions()
383 : @*/
384 517 : PetscErrorCode EPSSetDimensions(EPS eps,PetscInt nev,PetscInt ncv,PetscInt mpd)
385 : {
386 517 : PetscFunctionBegin;
387 517 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
388 2068 : PetscValidLogicalCollectiveInt(eps,nev,2);
389 2068 : PetscValidLogicalCollectiveInt(eps,ncv,3);
390 2068 : PetscValidLogicalCollectiveInt(eps,mpd,4);
391 517 : PetscCheck(nev>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
392 517 : eps->nev = nev;
393 517 : if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
394 364 : eps->ncv = PETSC_DEFAULT;
395 : } else {
396 153 : PetscCheck(ncv>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
397 153 : eps->ncv = ncv;
398 : }
399 517 : if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
400 466 : eps->mpd = PETSC_DEFAULT;
401 : } else {
402 51 : PetscCheck(mpd>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
403 51 : eps->mpd = mpd;
404 : }
405 517 : eps->state = EPS_STATE_INITIAL;
406 517 : PetscFunctionReturn(PETSC_SUCCESS);
407 : }
408 :
409 : /*@
410 : EPSSetWhichEigenpairs - Specifies which portion of the spectrum is
411 : to be sought.
412 :
413 : Logically Collective
414 :
415 : Input Parameters:
416 : + eps - eigensolver context obtained from EPSCreate()
417 : - which - the portion of the spectrum to be sought
418 :
419 : Options Database Keys:
420 : + -eps_largest_magnitude - Sets largest eigenvalues in magnitude
421 : . -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
422 : . -eps_largest_real - Sets largest real parts
423 : . -eps_smallest_real - Sets smallest real parts
424 : . -eps_largest_imaginary - Sets largest imaginary parts
425 : . -eps_smallest_imaginary - Sets smallest imaginary parts
426 : . -eps_target_magnitude - Sets eigenvalues closest to target
427 : . -eps_target_real - Sets real parts closest to target
428 : . -eps_target_imaginary - Sets imaginary parts closest to target
429 : - -eps_all - Sets all eigenvalues in an interval or region
430 :
431 : Notes:
432 : The parameter 'which' can have one of these values
433 :
434 : + EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
435 : . EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
436 : . EPS_LARGEST_REAL - largest real parts
437 : . EPS_SMALLEST_REAL - smallest real parts
438 : . EPS_LARGEST_IMAGINARY - largest imaginary parts
439 : . EPS_SMALLEST_IMAGINARY - smallest imaginary parts
440 : . EPS_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
441 : . EPS_TARGET_REAL - eigenvalues with real part closest to target
442 : . EPS_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
443 : . EPS_ALL - all eigenvalues contained in a given interval or region
444 : - EPS_WHICH_USER - user defined ordering set with EPSSetEigenvalueComparison()
445 :
446 : Not all eigensolvers implemented in EPS account for all the possible values
447 : stated above. Also, some values make sense only for certain types of
448 : problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY
449 : and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part
450 : for eigenvalue selection.
451 :
452 : The target is a scalar value provided with EPSSetTarget().
453 :
454 : The criterion EPS_TARGET_IMAGINARY is available only in case PETSc and
455 : SLEPc have been built with complex scalars.
456 :
457 : EPS_ALL is intended for use in combination with an interval (see
458 : EPSSetInterval()), when all eigenvalues within the interval are requested,
459 : or in the context of the CISS solver for computing all eigenvalues in a region.
460 : In those cases, the number of eigenvalues is unknown, so the nev parameter
461 : has a different sense, see EPSSetDimensions().
462 :
463 : Level: intermediate
464 :
465 : .seealso: EPSGetWhichEigenpairs(), EPSSetTarget(), EPSSetInterval(),
466 : EPSSetDimensions(), EPSSetEigenvalueComparison(), EPSWhich
467 : @*/
468 580 : PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)
469 : {
470 580 : PetscFunctionBegin;
471 580 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
472 2320 : PetscValidLogicalCollectiveEnum(eps,which,2);
473 580 : switch (which) {
474 580 : case EPS_LARGEST_MAGNITUDE:
475 : case EPS_SMALLEST_MAGNITUDE:
476 : case EPS_LARGEST_REAL:
477 : case EPS_SMALLEST_REAL:
478 : case EPS_LARGEST_IMAGINARY:
479 : case EPS_SMALLEST_IMAGINARY:
480 : case EPS_TARGET_MAGNITUDE:
481 : case EPS_TARGET_REAL:
482 : #if defined(PETSC_USE_COMPLEX)
483 : case EPS_TARGET_IMAGINARY:
484 : #endif
485 : case EPS_ALL:
486 : case EPS_WHICH_USER:
487 580 : if (eps->which != which) {
488 502 : eps->state = EPS_STATE_INITIAL;
489 502 : eps->which = which;
490 : }
491 580 : break;
492 : #if !defined(PETSC_USE_COMPLEX)
493 : case EPS_TARGET_IMAGINARY:
494 : SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"EPS_TARGET_IMAGINARY can be used only with complex scalars");
495 : #endif
496 0 : default:
497 0 : SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
498 : }
499 580 : PetscFunctionReturn(PETSC_SUCCESS);
500 : }
501 :
502 : /*@
503 : EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be
504 : sought.
505 :
506 : Not Collective
507 :
508 : Input Parameter:
509 : . eps - eigensolver context obtained from EPSCreate()
510 :
511 : Output Parameter:
512 : . which - the portion of the spectrum to be sought
513 :
514 : Notes:
515 : See EPSSetWhichEigenpairs() for possible values of 'which'.
516 :
517 : Level: intermediate
518 :
519 : .seealso: EPSSetWhichEigenpairs(), EPSWhich
520 : @*/
521 1 : PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)
522 : {
523 1 : PetscFunctionBegin;
524 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
525 1 : PetscAssertPointer(which,2);
526 1 : *which = eps->which;
527 1 : PetscFunctionReturn(PETSC_SUCCESS);
528 : }
529 :
530 : /*@C
531 : EPSSetEigenvalueComparison - Specifies the eigenvalue comparison function
532 : when EPSSetWhichEigenpairs() is set to EPS_WHICH_USER.
533 :
534 : Logically Collective
535 :
536 : Input Parameters:
537 : + eps - eigensolver context obtained from EPSCreate()
538 : . func - a pointer to the comparison function
539 : - ctx - a context pointer (the last parameter to the comparison function)
540 :
541 : Calling sequence of func:
542 : $ PetscErrorCode func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
543 : + ar - real part of the 1st eigenvalue
544 : . ai - imaginary part of the 1st eigenvalue
545 : . br - real part of the 2nd eigenvalue
546 : . bi - imaginary part of the 2nd eigenvalue
547 : . res - result of comparison
548 : - ctx - optional context, as set by EPSSetEigenvalueComparison()
549 :
550 : Note:
551 : The returning parameter 'res' can be
552 : + negative - if the 1st eigenvalue is preferred to the 2st one
553 : . zero - if both eigenvalues are equally preferred
554 : - positive - if the 2st eigenvalue is preferred to the 1st one
555 :
556 : Level: advanced
557 :
558 : .seealso: EPSSetWhichEigenpairs(), EPSWhich
559 : @*/
560 38 : PetscErrorCode EPSSetEigenvalueComparison(EPS eps,PetscErrorCode (*func)(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx),void* ctx)
561 : {
562 38 : PetscFunctionBegin;
563 38 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
564 38 : eps->sc->comparison = func;
565 38 : eps->sc->comparisonctx = ctx;
566 38 : eps->which = EPS_WHICH_USER;
567 38 : PetscFunctionReturn(PETSC_SUCCESS);
568 : }
569 :
570 : /*@C
571 : EPSSetArbitrarySelection - Specifies a function intended to look for
572 : eigenvalues according to an arbitrary selection criterion. This criterion
573 : can be based on a computation involving the current eigenvector approximation.
574 :
575 : Logically Collective
576 :
577 : Input Parameters:
578 : + eps - eigensolver context obtained from EPSCreate()
579 : . func - a pointer to the evaluation function
580 : - ctx - a context pointer (the last parameter to the evaluation function)
581 :
582 : Calling sequence of func:
583 : $ PetscErrorCode func(PetscScalar er,PetscScalar ei,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx)
584 : + er - real part of the current eigenvalue approximation
585 : . ei - imaginary part of the current eigenvalue approximation
586 : . xr - real part of the current eigenvector approximation
587 : . xi - imaginary part of the current eigenvector approximation
588 : . rr - result of evaluation (real part)
589 : . ri - result of evaluation (imaginary part)
590 : - ctx - optional context, as set by EPSSetArbitrarySelection()
591 :
592 : Notes:
593 : This provides a mechanism to select eigenpairs by evaluating a user-defined
594 : function. When a function has been provided, the default selection based on
595 : sorting the eigenvalues is replaced by the sorting of the results of this
596 : function (with the same sorting criterion given in EPSSetWhichEigenpairs()).
597 :
598 : For instance, suppose you want to compute those eigenvectors that maximize
599 : a certain computable expression. Then implement the computation using
600 : the arguments xr and xi, and return the result in rr. Then set the standard
601 : sorting by magnitude so that the eigenpair with largest value of rr is
602 : selected.
603 :
604 : This evaluation function is collective, that is, all processes call it and
605 : it can use collective operations; furthermore, the computed result must
606 : be the same in all processes.
607 :
608 : The result of func is expressed as a complex number so that it is possible to
609 : use the standard eigenvalue sorting functions, but normally only rr is used.
610 : Set ri to zero unless it is meaningful in your application.
611 :
612 : Level: advanced
613 :
614 : .seealso: EPSSetWhichEigenpairs()
615 : @*/
616 8 : PetscErrorCode EPSSetArbitrarySelection(EPS eps,PetscErrorCode (*func)(PetscScalar er,PetscScalar ei,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx),void* ctx)
617 : {
618 8 : PetscFunctionBegin;
619 8 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
620 8 : eps->arbitrary = func;
621 8 : eps->arbitraryctx = ctx;
622 8 : eps->state = EPS_STATE_INITIAL;
623 8 : PetscFunctionReturn(PETSC_SUCCESS);
624 : }
625 :
626 : /*@C
627 : EPSSetConvergenceTestFunction - Sets a function to compute the error estimate
628 : used in the convergence test.
629 :
630 : Logically Collective
631 :
632 : Input Parameters:
633 : + eps - eigensolver context obtained from EPSCreate()
634 : . func - a pointer to the convergence test function
635 : . ctx - context for private data for the convergence routine (may be null)
636 : - destroy - a routine for destroying the context (may be null)
637 :
638 : Calling sequence of func:
639 : $ PetscErrorCode func(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
640 : + eps - eigensolver context obtained from EPSCreate()
641 : . eigr - real part of the eigenvalue
642 : . eigi - imaginary part of the eigenvalue
643 : . res - residual norm associated to the eigenpair
644 : . errest - (output) computed error estimate
645 : - ctx - optional context, as set by EPSSetConvergenceTestFunction()
646 :
647 : Note:
648 : If the error estimate returned by the convergence test function is less than
649 : the tolerance, then the eigenvalue is accepted as converged.
650 :
651 : Level: advanced
652 :
653 : .seealso: EPSSetConvergenceTest(), EPSSetTolerances()
654 : @*/
655 32 : PetscErrorCode EPSSetConvergenceTestFunction(EPS eps,PetscErrorCode (*func)(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx),void* ctx,PetscErrorCode (*destroy)(void*))
656 : {
657 32 : PetscFunctionBegin;
658 32 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
659 32 : if (eps->convergeddestroy) PetscCall((*eps->convergeddestroy)(eps->convergedctx));
660 32 : eps->convergeduser = func;
661 32 : eps->convergeddestroy = destroy;
662 32 : eps->convergedctx = ctx;
663 32 : if (func == EPSConvergedRelative) eps->conv = EPS_CONV_REL;
664 32 : else if (func == EPSConvergedNorm) eps->conv = EPS_CONV_NORM;
665 32 : else if (func == EPSConvergedAbsolute) eps->conv = EPS_CONV_ABS;
666 : else {
667 32 : eps->conv = EPS_CONV_USER;
668 32 : eps->converged = eps->convergeduser;
669 : }
670 32 : PetscFunctionReturn(PETSC_SUCCESS);
671 : }
672 :
673 : /*@
674 : EPSSetConvergenceTest - Specifies how to compute the error estimate
675 : used in the convergence test.
676 :
677 : Logically Collective
678 :
679 : Input Parameters:
680 : + eps - eigensolver context obtained from EPSCreate()
681 : - conv - the type of convergence test
682 :
683 : Options Database Keys:
684 : + -eps_conv_abs - Sets the absolute convergence test
685 : . -eps_conv_rel - Sets the convergence test relative to the eigenvalue
686 : . -eps_conv_norm - Sets the convergence test relative to the matrix norms
687 : - -eps_conv_user - Selects the user-defined convergence test
688 :
689 : Note:
690 : The parameter 'conv' can have one of these values
691 : + EPS_CONV_ABS - absolute error ||r||
692 : . EPS_CONV_REL - error relative to the eigenvalue l, ||r||/|l|
693 : . EPS_CONV_NORM - error relative to the matrix norms, ||r||/(||A||+|l|*||B||)
694 : - EPS_CONV_USER - function set by EPSSetConvergenceTestFunction()
695 :
696 : Level: intermediate
697 :
698 : .seealso: EPSGetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSSetStoppingTest(), EPSConv
699 : @*/
700 153 : PetscErrorCode EPSSetConvergenceTest(EPS eps,EPSConv conv)
701 : {
702 153 : PetscFunctionBegin;
703 153 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
704 612 : PetscValidLogicalCollectiveEnum(eps,conv,2);
705 153 : switch (conv) {
706 19 : case EPS_CONV_ABS: eps->converged = EPSConvergedAbsolute; break;
707 84 : case EPS_CONV_REL: eps->converged = EPSConvergedRelative; break;
708 50 : case EPS_CONV_NORM: eps->converged = EPSConvergedNorm; break;
709 0 : case EPS_CONV_USER:
710 0 : PetscCheck(eps->convergeduser,PetscObjectComm((PetscObject)eps),PETSC_ERR_ORDER,"Must call EPSSetConvergenceTestFunction() first");
711 0 : eps->converged = eps->convergeduser;
712 0 : break;
713 0 : default:
714 0 : SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
715 : }
716 153 : eps->conv = conv;
717 153 : PetscFunctionReturn(PETSC_SUCCESS);
718 : }
719 :
720 : /*@
721 : EPSGetConvergenceTest - Gets the method used to compute the error estimate
722 : used in the convergence test.
723 :
724 : Not Collective
725 :
726 : Input Parameters:
727 : . eps - eigensolver context obtained from EPSCreate()
728 :
729 : Output Parameters:
730 : . conv - the type of convergence test
731 :
732 : Level: intermediate
733 :
734 : .seealso: EPSSetConvergenceTest(), EPSConv
735 : @*/
736 1 : PetscErrorCode EPSGetConvergenceTest(EPS eps,EPSConv *conv)
737 : {
738 1 : PetscFunctionBegin;
739 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
740 1 : PetscAssertPointer(conv,2);
741 1 : *conv = eps->conv;
742 1 : PetscFunctionReturn(PETSC_SUCCESS);
743 : }
744 :
745 : /*@C
746 : EPSSetStoppingTestFunction - Sets a function to decide when to stop the outer
747 : iteration of the eigensolver.
748 :
749 : Logically Collective
750 :
751 : Input Parameters:
752 : + eps - eigensolver context obtained from EPSCreate()
753 : . func - pointer to the stopping test function
754 : . ctx - context for private data for the stopping routine (may be null)
755 : - destroy - a routine for destroying the context (may be null)
756 :
757 : Calling sequence of func:
758 : $ PetscErrorCode func(EPS eps,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,EPSConvergedReason *reason,void *ctx)
759 : + eps - eigensolver context obtained from EPSCreate()
760 : . its - current number of iterations
761 : . max_it - maximum number of iterations
762 : . nconv - number of currently converged eigenpairs
763 : . nev - number of requested eigenpairs
764 : . reason - (output) result of the stopping test
765 : - ctx - optional context, as set by EPSSetStoppingTestFunction()
766 :
767 : Note:
768 : Normal usage is to first call the default routine EPSStoppingBasic() and then
769 : set reason to EPS_CONVERGED_USER if some user-defined conditions have been
770 : met. To let the eigensolver continue iterating, the result must be left as
771 : EPS_CONVERGED_ITERATING.
772 :
773 : Level: advanced
774 :
775 : .seealso: EPSSetStoppingTest(), EPSStoppingBasic()
776 : @*/
777 2 : PetscErrorCode EPSSetStoppingTestFunction(EPS eps,PetscErrorCode (*func)(EPS eps,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,EPSConvergedReason *reason,void *ctx),void* ctx,PetscErrorCode (*destroy)(void*))
778 : {
779 2 : PetscFunctionBegin;
780 2 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
781 2 : if (eps->stoppingdestroy) PetscCall((*eps->stoppingdestroy)(eps->stoppingctx));
782 2 : eps->stoppinguser = func;
783 2 : eps->stoppingdestroy = destroy;
784 2 : eps->stoppingctx = ctx;
785 2 : if (func == EPSStoppingBasic) eps->stop = EPS_STOP_BASIC;
786 : else {
787 2 : eps->stop = EPS_STOP_USER;
788 2 : eps->stopping = eps->stoppinguser;
789 : }
790 2 : PetscFunctionReturn(PETSC_SUCCESS);
791 : }
792 :
793 : /*@
794 : EPSSetStoppingTest - Specifies how to decide the termination of the outer
795 : loop of the eigensolver.
796 :
797 : Logically Collective
798 :
799 : Input Parameters:
800 : + eps - eigensolver context obtained from EPSCreate()
801 : - stop - the type of stopping test
802 :
803 : Options Database Keys:
804 : + -eps_stop_basic - Sets the default stopping test
805 : - -eps_stop_user - Selects the user-defined stopping test
806 :
807 : Note:
808 : The parameter 'stop' can have one of these values
809 : + EPS_STOP_BASIC - default stopping test
810 : - EPS_STOP_USER - function set by EPSSetStoppingTestFunction()
811 :
812 : Level: advanced
813 :
814 : .seealso: EPSGetStoppingTest(), EPSSetStoppingTestFunction(), EPSSetConvergenceTest(), EPSStop
815 : @*/
816 1 : PetscErrorCode EPSSetStoppingTest(EPS eps,EPSStop stop)
817 : {
818 1 : PetscFunctionBegin;
819 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
820 4 : PetscValidLogicalCollectiveEnum(eps,stop,2);
821 1 : switch (stop) {
822 1 : case EPS_STOP_BASIC: eps->stopping = EPSStoppingBasic; break;
823 0 : case EPS_STOP_USER:
824 0 : PetscCheck(eps->stoppinguser,PetscObjectComm((PetscObject)eps),PETSC_ERR_ORDER,"Must call EPSSetStoppingTestFunction() first");
825 0 : eps->stopping = eps->stoppinguser;
826 0 : break;
827 0 : default:
828 0 : SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
829 : }
830 1 : eps->stop = stop;
831 1 : PetscFunctionReturn(PETSC_SUCCESS);
832 : }
833 :
834 : /*@
835 : EPSGetStoppingTest - Gets the method used to decide the termination of the outer
836 : loop of the eigensolver.
837 :
838 : Not Collective
839 :
840 : Input Parameters:
841 : . eps - eigensolver context obtained from EPSCreate()
842 :
843 : Output Parameters:
844 : . stop - the type of stopping test
845 :
846 : Level: advanced
847 :
848 : .seealso: EPSSetStoppingTest(), EPSStop
849 : @*/
850 1 : PetscErrorCode EPSGetStoppingTest(EPS eps,EPSStop *stop)
851 : {
852 1 : PetscFunctionBegin;
853 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
854 1 : PetscAssertPointer(stop,2);
855 1 : *stop = eps->stop;
856 1 : PetscFunctionReturn(PETSC_SUCCESS);
857 : }
858 :
859 : /*@
860 : EPSSetProblemType - Specifies the type of the eigenvalue problem.
861 :
862 : Logically Collective
863 :
864 : Input Parameters:
865 : + eps - the eigensolver context
866 : - type - a known type of eigenvalue problem
867 :
868 : Options Database Keys:
869 : + -eps_hermitian - Hermitian eigenvalue problem
870 : . -eps_gen_hermitian - generalized Hermitian eigenvalue problem
871 : . -eps_non_hermitian - non-Hermitian eigenvalue problem
872 : . -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
873 : . -eps_pos_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
874 : with positive semi-definite B
875 : - -eps_gen_indefinite - generalized Hermitian-indefinite eigenvalue problem
876 :
877 : Notes:
878 : Allowed values for the problem type are Hermitian (EPS_HEP), non-Hermitian
879 : (EPS_NHEP), generalized Hermitian (EPS_GHEP), generalized non-Hermitian
880 : (EPS_GNHEP), generalized non-Hermitian with positive semi-definite B
881 : (EPS_PGNHEP), and generalized Hermitian-indefinite (EPS_GHIEP).
882 :
883 : This function must be used to instruct SLEPc to exploit symmetry. If no
884 : problem type is specified, by default a non-Hermitian problem is assumed
885 : (either standard or generalized). If the user knows that the problem is
886 : Hermitian (i.e. A=A^H) or generalized Hermitian (i.e. A=A^H, B=B^H, and
887 : B positive definite) then it is recommended to set the problem type so
888 : that eigensolver can exploit these properties.
889 :
890 : Level: intermediate
891 :
892 : .seealso: EPSSetOperators(), EPSSetType(), EPSGetProblemType(), EPSProblemType
893 : @*/
894 643 : PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)
895 : {
896 643 : PetscFunctionBegin;
897 643 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
898 2572 : PetscValidLogicalCollectiveEnum(eps,type,2);
899 643 : if (type == eps->problem_type) PetscFunctionReturn(PETSC_SUCCESS);
900 618 : switch (type) {
901 248 : case EPS_HEP:
902 248 : eps->isgeneralized = PETSC_FALSE;
903 248 : eps->ishermitian = PETSC_TRUE;
904 248 : eps->ispositive = PETSC_FALSE;
905 248 : break;
906 168 : case EPS_NHEP:
907 168 : eps->isgeneralized = PETSC_FALSE;
908 168 : eps->ishermitian = PETSC_FALSE;
909 168 : eps->ispositive = PETSC_FALSE;
910 168 : break;
911 125 : case EPS_GHEP:
912 125 : eps->isgeneralized = PETSC_TRUE;
913 125 : eps->ishermitian = PETSC_TRUE;
914 125 : eps->ispositive = PETSC_TRUE;
915 125 : break;
916 59 : case EPS_GNHEP:
917 59 : eps->isgeneralized = PETSC_TRUE;
918 59 : eps->ishermitian = PETSC_FALSE;
919 59 : eps->ispositive = PETSC_FALSE;
920 59 : break;
921 2 : case EPS_PGNHEP:
922 2 : eps->isgeneralized = PETSC_TRUE;
923 2 : eps->ishermitian = PETSC_FALSE;
924 2 : eps->ispositive = PETSC_TRUE;
925 2 : break;
926 16 : case EPS_GHIEP:
927 16 : eps->isgeneralized = PETSC_TRUE;
928 16 : eps->ishermitian = PETSC_TRUE;
929 16 : eps->ispositive = PETSC_FALSE;
930 16 : break;
931 0 : default:
932 0 : SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
933 : }
934 618 : eps->problem_type = type;
935 618 : eps->state = EPS_STATE_INITIAL;
936 618 : PetscFunctionReturn(PETSC_SUCCESS);
937 : }
938 :
939 : /*@
940 : EPSGetProblemType - Gets the problem type from the EPS object.
941 :
942 : Not Collective
943 :
944 : Input Parameter:
945 : . eps - the eigensolver context
946 :
947 : Output Parameter:
948 : . type - the problem type
949 :
950 : Level: intermediate
951 :
952 : .seealso: EPSSetProblemType(), EPSProblemType
953 : @*/
954 59 : PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)
955 : {
956 59 : PetscFunctionBegin;
957 59 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
958 59 : PetscAssertPointer(type,2);
959 59 : *type = eps->problem_type;
960 59 : PetscFunctionReturn(PETSC_SUCCESS);
961 : }
962 :
963 : /*@
964 : EPSSetExtraction - Specifies the type of extraction technique to be employed
965 : by the eigensolver.
966 :
967 : Logically Collective
968 :
969 : Input Parameters:
970 : + eps - the eigensolver context
971 : - extr - a known type of extraction
972 :
973 : Options Database Keys:
974 : + -eps_ritz - Rayleigh-Ritz extraction
975 : . -eps_harmonic - harmonic Ritz extraction
976 : . -eps_harmonic_relative - harmonic Ritz extraction relative to the eigenvalue
977 : . -eps_harmonic_right - harmonic Ritz extraction for rightmost eigenvalues
978 : . -eps_harmonic_largest - harmonic Ritz extraction for largest magnitude
979 : (without target)
980 : . -eps_refined - refined Ritz extraction
981 : - -eps_refined_harmonic - refined harmonic Ritz extraction
982 :
983 : Notes:
984 : Not all eigensolvers support all types of extraction. See the SLEPc
985 : Users Manual for details.
986 :
987 : By default, a standard Rayleigh-Ritz extraction is used. Other extractions
988 : may be useful when computing interior eigenvalues.
989 :
990 : Harmonic-type extractions are used in combination with a 'target'.
991 :
992 : Level: advanced
993 :
994 : .seealso: EPSSetTarget(), EPSGetExtraction(), EPSExtraction
995 : @*/
996 14 : PetscErrorCode EPSSetExtraction(EPS eps,EPSExtraction extr)
997 : {
998 14 : PetscFunctionBegin;
999 14 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1000 56 : PetscValidLogicalCollectiveEnum(eps,extr,2);
1001 14 : eps->extraction = extr;
1002 14 : PetscFunctionReturn(PETSC_SUCCESS);
1003 : }
1004 :
1005 : /*@
1006 : EPSGetExtraction - Gets the extraction type used by the EPS object.
1007 :
1008 : Not Collective
1009 :
1010 : Input Parameter:
1011 : . eps - the eigensolver context
1012 :
1013 : Output Parameter:
1014 : . extr - name of extraction type
1015 :
1016 : Level: advanced
1017 :
1018 : .seealso: EPSSetExtraction(), EPSExtraction
1019 : @*/
1020 2 : PetscErrorCode EPSGetExtraction(EPS eps,EPSExtraction *extr)
1021 : {
1022 2 : PetscFunctionBegin;
1023 2 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1024 2 : PetscAssertPointer(extr,2);
1025 2 : *extr = eps->extraction;
1026 2 : PetscFunctionReturn(PETSC_SUCCESS);
1027 : }
1028 :
1029 : /*@
1030 : EPSSetBalance - Specifies the balancing technique to be employed by the
1031 : eigensolver, and some parameters associated to it.
1032 :
1033 : Logically Collective
1034 :
1035 : Input Parameters:
1036 : + eps - the eigensolver context
1037 : . bal - the balancing method, one of EPS_BALANCE_NONE, EPS_BALANCE_ONESIDE,
1038 : EPS_BALANCE_TWOSIDE, or EPS_BALANCE_USER
1039 : . its - number of iterations of the balancing algorithm
1040 : - cutoff - cutoff value
1041 :
1042 : Options Database Keys:
1043 : + -eps_balance <method> - the balancing method, where <method> is one of
1044 : 'none', 'oneside', 'twoside', or 'user'
1045 : . -eps_balance_its <its> - number of iterations
1046 : - -eps_balance_cutoff <cutoff> - cutoff value
1047 :
1048 : Notes:
1049 : When balancing is enabled, the solver works implicitly with matrix DAD^-1,
1050 : where D is an appropriate diagonal matrix. This improves the accuracy of
1051 : the computed results in some cases. See the SLEPc Users Manual for details.
1052 :
1053 : Balancing makes sense only for non-Hermitian problems when the required
1054 : precision is high (i.e. a small tolerance such as 1e-15).
1055 :
1056 : By default, balancing is disabled. The two-sided method is much more
1057 : effective than the one-sided counterpart, but it requires the system
1058 : matrices to have the MatMultTranspose operation defined.
1059 :
1060 : The parameter 'its' is the number of iterations performed by the method. The
1061 : cutoff value is used only in the two-side variant. Use PETSC_DEFAULT to assign
1062 : a reasonably good value.
1063 :
1064 : User-defined balancing is allowed provided that the corresponding matrix
1065 : is set via STSetBalanceMatrix.
1066 :
1067 : Level: intermediate
1068 :
1069 : .seealso: EPSGetBalance(), EPSBalance, STSetBalanceMatrix()
1070 : @*/
1071 16 : PetscErrorCode EPSSetBalance(EPS eps,EPSBalance bal,PetscInt its,PetscReal cutoff)
1072 : {
1073 16 : PetscFunctionBegin;
1074 16 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1075 64 : PetscValidLogicalCollectiveEnum(eps,bal,2);
1076 64 : PetscValidLogicalCollectiveInt(eps,its,3);
1077 64 : PetscValidLogicalCollectiveReal(eps,cutoff,4);
1078 16 : switch (bal) {
1079 16 : case EPS_BALANCE_NONE:
1080 : case EPS_BALANCE_ONESIDE:
1081 : case EPS_BALANCE_TWOSIDE:
1082 : case EPS_BALANCE_USER:
1083 16 : if (eps->balance != bal) {
1084 14 : eps->state = EPS_STATE_INITIAL;
1085 14 : eps->balance = bal;
1086 : }
1087 16 : break;
1088 0 : default:
1089 0 : SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid value of argument 'bal'");
1090 : }
1091 16 : if (its==PETSC_DECIDE || its==PETSC_DEFAULT) eps->balance_its = 5;
1092 : else {
1093 16 : PetscCheck(its>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be > 0");
1094 16 : eps->balance_its = its;
1095 : }
1096 16 : if (cutoff==(PetscReal)PETSC_DECIDE || cutoff==(PetscReal)PETSC_DEFAULT) eps->balance_cutoff = 1e-8;
1097 : else {
1098 16 : PetscCheck(cutoff>0.0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of cutoff. Must be > 0");
1099 16 : eps->balance_cutoff = cutoff;
1100 : }
1101 16 : PetscFunctionReturn(PETSC_SUCCESS);
1102 : }
1103 :
1104 : /*@C
1105 : EPSGetBalance - Gets the balancing type used by the EPS object, and the
1106 : associated parameters.
1107 :
1108 : Not Collective
1109 :
1110 : Input Parameter:
1111 : . eps - the eigensolver context
1112 :
1113 : Output Parameters:
1114 : + bal - the balancing method
1115 : . its - number of iterations of the balancing algorithm
1116 : - cutoff - cutoff value
1117 :
1118 : Level: intermediate
1119 :
1120 : Note:
1121 : The user can specify NULL for any parameter that is not needed.
1122 :
1123 : .seealso: EPSSetBalance(), EPSBalance
1124 : @*/
1125 1 : PetscErrorCode EPSGetBalance(EPS eps,EPSBalance *bal,PetscInt *its,PetscReal *cutoff)
1126 : {
1127 1 : PetscFunctionBegin;
1128 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1129 1 : if (bal) *bal = eps->balance;
1130 1 : if (its) *its = eps->balance_its;
1131 1 : if (cutoff) *cutoff = eps->balance_cutoff;
1132 1 : PetscFunctionReturn(PETSC_SUCCESS);
1133 : }
1134 :
1135 : /*@
1136 : EPSSetTwoSided - Sets the solver to use a two-sided variant so that left
1137 : eigenvectors are also computed.
1138 :
1139 : Logically Collective
1140 :
1141 : Input Parameters:
1142 : + eps - the eigensolver context
1143 : - twosided - whether the two-sided variant is to be used or not
1144 :
1145 : Options Database Keys:
1146 : . -eps_two_sided <boolean> - Sets/resets the twosided flag
1147 :
1148 : Notes:
1149 : If the user sets twosided=PETSC_TRUE then the solver uses a variant of
1150 : the algorithm that computes both right and left eigenvectors. This is
1151 : usually much more costly. This option is not available in all solvers.
1152 :
1153 : When using two-sided solvers, the problem matrices must have both the
1154 : MatMult and MatMultTranspose operations defined.
1155 :
1156 : Level: advanced
1157 :
1158 : .seealso: EPSGetTwoSided(), EPSGetLeftEigenvector()
1159 : @*/
1160 26 : PetscErrorCode EPSSetTwoSided(EPS eps,PetscBool twosided)
1161 : {
1162 26 : PetscFunctionBegin;
1163 26 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1164 104 : PetscValidLogicalCollectiveBool(eps,twosided,2);
1165 26 : if (twosided!=eps->twosided) {
1166 17 : eps->twosided = twosided;
1167 17 : eps->state = EPS_STATE_INITIAL;
1168 : }
1169 26 : PetscFunctionReturn(PETSC_SUCCESS);
1170 : }
1171 :
1172 : /*@
1173 : EPSGetTwoSided - Returns the flag indicating whether a two-sided variant
1174 : of the algorithm is being used or not.
1175 :
1176 : Not Collective
1177 :
1178 : Input Parameter:
1179 : . eps - the eigensolver context
1180 :
1181 : Output Parameter:
1182 : . twosided - the returned flag
1183 :
1184 : Level: advanced
1185 :
1186 : .seealso: EPSSetTwoSided()
1187 : @*/
1188 5 : PetscErrorCode EPSGetTwoSided(EPS eps,PetscBool *twosided)
1189 : {
1190 5 : PetscFunctionBegin;
1191 5 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1192 5 : PetscAssertPointer(twosided,2);
1193 5 : *twosided = eps->twosided;
1194 5 : PetscFunctionReturn(PETSC_SUCCESS);
1195 : }
1196 :
1197 : /*@
1198 : EPSSetTrueResidual - Specifies if the solver must compute the true residual
1199 : explicitly or not.
1200 :
1201 : Logically Collective
1202 :
1203 : Input Parameters:
1204 : + eps - the eigensolver context
1205 : - trueres - whether true residuals are required or not
1206 :
1207 : Options Database Keys:
1208 : . -eps_true_residual <boolean> - Sets/resets the boolean flag 'trueres'
1209 :
1210 : Notes:
1211 : If the user sets trueres=PETSC_TRUE then the solver explicitly computes
1212 : the true residual for each eigenpair approximation, and uses it for
1213 : convergence testing. Computing the residual is usually an expensive
1214 : operation. Some solvers (e.g., Krylov solvers) can avoid this computation
1215 : by using a cheap estimate of the residual norm, but this may sometimes
1216 : give inaccurate results (especially if a spectral transform is being
1217 : used). On the contrary, preconditioned eigensolvers (e.g., Davidson solvers)
1218 : do rely on computing the true residual, so this option is irrelevant for them.
1219 :
1220 : Level: advanced
1221 :
1222 : .seealso: EPSGetTrueResidual()
1223 : @*/
1224 9 : PetscErrorCode EPSSetTrueResidual(EPS eps,PetscBool trueres)
1225 : {
1226 9 : PetscFunctionBegin;
1227 9 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1228 36 : PetscValidLogicalCollectiveBool(eps,trueres,2);
1229 9 : eps->trueres = trueres;
1230 9 : PetscFunctionReturn(PETSC_SUCCESS);
1231 : }
1232 :
1233 : /*@
1234 : EPSGetTrueResidual - Returns the flag indicating whether true
1235 : residuals must be computed explicitly or not.
1236 :
1237 : Not Collective
1238 :
1239 : Input Parameter:
1240 : . eps - the eigensolver context
1241 :
1242 : Output Parameter:
1243 : . trueres - the returned flag
1244 :
1245 : Level: advanced
1246 :
1247 : .seealso: EPSSetTrueResidual()
1248 : @*/
1249 8 : PetscErrorCode EPSGetTrueResidual(EPS eps,PetscBool *trueres)
1250 : {
1251 8 : PetscFunctionBegin;
1252 8 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1253 8 : PetscAssertPointer(trueres,2);
1254 8 : *trueres = eps->trueres;
1255 8 : PetscFunctionReturn(PETSC_SUCCESS);
1256 : }
1257 :
1258 : /*@
1259 : EPSSetTrackAll - Specifies if the solver must compute the residual norm of all
1260 : approximate eigenpairs or not.
1261 :
1262 : Logically Collective
1263 :
1264 : Input Parameters:
1265 : + eps - the eigensolver context
1266 : - trackall - whether to compute all residuals or not
1267 :
1268 : Notes:
1269 : If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
1270 : the residual norm for each eigenpair approximation. Computing the residual is
1271 : usually an expensive operation and solvers commonly compute only the residual
1272 : associated to the first unconverged eigenpair.
1273 :
1274 : The option '-eps_monitor_all' automatically activates this option.
1275 :
1276 : Level: developer
1277 :
1278 : .seealso: EPSGetTrackAll()
1279 : @*/
1280 134 : PetscErrorCode EPSSetTrackAll(EPS eps,PetscBool trackall)
1281 : {
1282 134 : PetscFunctionBegin;
1283 134 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1284 536 : PetscValidLogicalCollectiveBool(eps,trackall,2);
1285 134 : eps->trackall = trackall;
1286 134 : PetscFunctionReturn(PETSC_SUCCESS);
1287 : }
1288 :
1289 : /*@
1290 : EPSGetTrackAll - Returns the flag indicating whether all residual norms must
1291 : be computed or not.
1292 :
1293 : Not Collective
1294 :
1295 : Input Parameter:
1296 : . eps - the eigensolver context
1297 :
1298 : Output Parameter:
1299 : . trackall - the returned flag
1300 :
1301 : Level: developer
1302 :
1303 : .seealso: EPSSetTrackAll()
1304 : @*/
1305 1 : PetscErrorCode EPSGetTrackAll(EPS eps,PetscBool *trackall)
1306 : {
1307 1 : PetscFunctionBegin;
1308 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1309 1 : PetscAssertPointer(trackall,2);
1310 1 : *trackall = eps->trackall;
1311 1 : PetscFunctionReturn(PETSC_SUCCESS);
1312 : }
1313 :
1314 : /*@
1315 : EPSSetPurify - Deactivate eigenvector purification (which is activated by default).
1316 :
1317 : Logically Collective
1318 :
1319 : Input Parameters:
1320 : + eps - the eigensolver context
1321 : - purify - whether purification is required or not
1322 :
1323 : Options Database Keys:
1324 : . -eps_purify <boolean> - Sets/resets the boolean flag 'purify'
1325 :
1326 : Notes:
1327 : By default, eigenvectors of generalized symmetric eigenproblems are purified
1328 : in order to purge directions in the nullspace of matrix B. If the user knows
1329 : that B is non-singular, then purification can be safely deactivated and some
1330 : computational cost is avoided (this is particularly important in interval computations).
1331 :
1332 : Level: intermediate
1333 :
1334 : .seealso: EPSGetPurify(), EPSSetInterval()
1335 : @*/
1336 4 : PetscErrorCode EPSSetPurify(EPS eps,PetscBool purify)
1337 : {
1338 4 : PetscFunctionBegin;
1339 4 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1340 16 : PetscValidLogicalCollectiveBool(eps,purify,2);
1341 4 : if (purify!=eps->purify) {
1342 2 : eps->purify = purify;
1343 2 : eps->state = EPS_STATE_INITIAL;
1344 : }
1345 4 : PetscFunctionReturn(PETSC_SUCCESS);
1346 : }
1347 :
1348 : /*@
1349 : EPSGetPurify - Returns the flag indicating whether purification is activated
1350 : or not.
1351 :
1352 : Not Collective
1353 :
1354 : Input Parameter:
1355 : . eps - the eigensolver context
1356 :
1357 : Output Parameter:
1358 : . purify - the returned flag
1359 :
1360 : Level: intermediate
1361 :
1362 : .seealso: EPSSetPurify()
1363 : @*/
1364 1 : PetscErrorCode EPSGetPurify(EPS eps,PetscBool *purify)
1365 : {
1366 1 : PetscFunctionBegin;
1367 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1368 1 : PetscAssertPointer(purify,2);
1369 1 : *purify = eps->purify;
1370 1 : PetscFunctionReturn(PETSC_SUCCESS);
1371 : }
1372 :
1373 : /*@C
1374 : EPSSetOptionsPrefix - Sets the prefix used for searching for all
1375 : EPS options in the database.
1376 :
1377 : Logically Collective
1378 :
1379 : Input Parameters:
1380 : + eps - the eigensolver context
1381 : - prefix - the prefix string to prepend to all EPS option requests
1382 :
1383 : Notes:
1384 : A hyphen (-) must NOT be given at the beginning of the prefix name.
1385 : The first character of all runtime options is AUTOMATICALLY the
1386 : hyphen.
1387 :
1388 : For example, to distinguish between the runtime options for two
1389 : different EPS contexts, one could call
1390 : .vb
1391 : EPSSetOptionsPrefix(eps1,"eig1_")
1392 : EPSSetOptionsPrefix(eps2,"eig2_")
1393 : .ve
1394 :
1395 : Level: advanced
1396 :
1397 : .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
1398 : @*/
1399 185 : PetscErrorCode EPSSetOptionsPrefix(EPS eps,const char *prefix)
1400 : {
1401 185 : PetscFunctionBegin;
1402 185 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1403 185 : if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
1404 185 : PetscCall(STSetOptionsPrefix(eps->st,prefix));
1405 185 : if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
1406 185 : PetscCall(BVSetOptionsPrefix(eps->V,prefix));
1407 185 : if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
1408 185 : PetscCall(DSSetOptionsPrefix(eps->ds,prefix));
1409 185 : if (!eps->rg) PetscCall(EPSGetRG(eps,&eps->rg));
1410 185 : PetscCall(RGSetOptionsPrefix(eps->rg,prefix));
1411 185 : PetscCall(PetscObjectSetOptionsPrefix((PetscObject)eps,prefix));
1412 185 : PetscFunctionReturn(PETSC_SUCCESS);
1413 : }
1414 :
1415 : /*@C
1416 : EPSAppendOptionsPrefix - Appends to the prefix used for searching for all
1417 : EPS options in the database.
1418 :
1419 : Logically Collective
1420 :
1421 : Input Parameters:
1422 : + eps - the eigensolver context
1423 : - prefix - the prefix string to prepend to all EPS option requests
1424 :
1425 : Notes:
1426 : A hyphen (-) must NOT be given at the beginning of the prefix name.
1427 : The first character of all runtime options is AUTOMATICALLY the hyphen.
1428 :
1429 : Level: advanced
1430 :
1431 : .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
1432 : @*/
1433 149 : PetscErrorCode EPSAppendOptionsPrefix(EPS eps,const char *prefix)
1434 : {
1435 149 : PetscFunctionBegin;
1436 149 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1437 149 : if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
1438 149 : PetscCall(STAppendOptionsPrefix(eps->st,prefix));
1439 149 : if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
1440 149 : PetscCall(BVAppendOptionsPrefix(eps->V,prefix));
1441 149 : if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
1442 149 : PetscCall(DSAppendOptionsPrefix(eps->ds,prefix));
1443 149 : if (!eps->rg) PetscCall(EPSGetRG(eps,&eps->rg));
1444 149 : PetscCall(RGAppendOptionsPrefix(eps->rg,prefix));
1445 149 : PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)eps,prefix));
1446 149 : PetscFunctionReturn(PETSC_SUCCESS);
1447 : }
1448 :
1449 : /*@C
1450 : EPSGetOptionsPrefix - Gets the prefix used for searching for all
1451 : EPS options in the database.
1452 :
1453 : Not Collective
1454 :
1455 : Input Parameters:
1456 : . eps - the eigensolver context
1457 :
1458 : Output Parameters:
1459 : . prefix - pointer to the prefix string used is returned
1460 :
1461 : Note:
1462 : On the Fortran side, the user should pass in a string 'prefix' of
1463 : sufficient length to hold the prefix.
1464 :
1465 : Level: advanced
1466 :
1467 : .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
1468 : @*/
1469 30 : PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])
1470 : {
1471 30 : PetscFunctionBegin;
1472 30 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
1473 30 : PetscAssertPointer(prefix,2);
1474 30 : PetscCall(PetscObjectGetOptionsPrefix((PetscObject)eps,prefix));
1475 30 : PetscFunctionReturn(PETSC_SUCCESS);
1476 : }
|