Actual source code: pepopts.c
slepc-3.21.0 2024-03-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: PEP routines related to options that can be set via the command-line
12: or procedurally
13: */
15: #include <slepc/private/pepimpl.h>
16: #include <petscdraw.h>
18: /*@C
19: PEPMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
20: indicated by the user.
22: Collective
24: Input Parameters:
25: + pep - the polynomial 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
31: Level: developer
33: .seealso: PEPMonitorSet(), PEPSetTrackAll()
34: @*/
35: PetscErrorCode PEPMonitorSetFromOptions(PEP pep,const char opt[],const char name[],void *ctx,PetscBool trackall)
36: {
37: PetscErrorCode (*mfunc)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
38: PetscErrorCode (*cfunc)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**);
39: PetscErrorCode (*dfunc)(PetscViewerAndFormat**);
40: PetscViewerAndFormat *vf;
41: PetscViewer viewer;
42: PetscViewerFormat format;
43: PetscViewerType vtype;
44: char key[PETSC_MAX_PATH_LEN];
45: PetscBool flg;
47: PetscFunctionBegin;
48: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,opt,&viewer,&format,&flg));
49: if (!flg) PetscFunctionReturn(PETSC_SUCCESS);
51: PetscCall(PetscViewerGetType(viewer,&vtype));
52: PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
53: PetscCall(PetscFunctionListFind(PEPMonitorList,key,&mfunc));
54: PetscCheck(mfunc,PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Specified viewer and format not supported");
55: PetscCall(PetscFunctionListFind(PEPMonitorCreateList,key,&cfunc));
56: PetscCall(PetscFunctionListFind(PEPMonitorDestroyList,key,&dfunc));
57: if (!cfunc) cfunc = PetscViewerAndFormatCreate_Internal;
58: if (!dfunc) dfunc = PetscViewerAndFormatDestroy;
60: PetscCall((*cfunc)(viewer,format,ctx,&vf));
61: PetscCall(PetscOptionsRestoreViewer(&viewer));
62: PetscCall(PEPMonitorSet(pep,mfunc,vf,(PetscErrorCode(*)(void **))dfunc));
63: if (trackall) PetscCall(PEPSetTrackAll(pep,PETSC_TRUE));
64: PetscFunctionReturn(PETSC_SUCCESS);
65: }
67: /*@
68: PEPSetFromOptions - Sets PEP options from the options database.
69: This routine must be called before PEPSetUp() if the user is to be
70: allowed to set the solver type.
72: Collective
74: Input Parameters:
75: . pep - the polynomial eigensolver context
77: Notes:
78: To see all options, run your program with the -help option.
80: Level: beginner
82: .seealso: PEPSetOptionsPrefix()
83: @*/
84: PetscErrorCode PEPSetFromOptions(PEP pep)
85: {
86: char type[256];
87: PetscBool set,flg,flg1,flg2,flg3,flg4,flg5;
88: PetscReal r,t,array[2]={0,0};
89: PetscScalar s;
90: PetscInt i,j,k;
91: PEPScale scale;
92: PEPRefine refine;
93: PEPRefineScheme scheme;
95: PetscFunctionBegin;
97: PetscCall(PEPRegisterAll());
98: PetscObjectOptionsBegin((PetscObject)pep);
99: PetscCall(PetscOptionsFList("-pep_type","Polynomial eigensolver method","PEPSetType",PEPList,(char*)(((PetscObject)pep)->type_name?((PetscObject)pep)->type_name:PEPTOAR),type,sizeof(type),&flg));
100: if (flg) PetscCall(PEPSetType(pep,type));
101: else if (!((PetscObject)pep)->type_name) PetscCall(PEPSetType(pep,PEPTOAR));
103: PetscCall(PetscOptionsBoolGroupBegin("-pep_general","General polynomial eigenvalue problem","PEPSetProblemType",&flg));
104: if (flg) PetscCall(PEPSetProblemType(pep,PEP_GENERAL));
105: PetscCall(PetscOptionsBoolGroup("-pep_hermitian","Hermitian polynomial eigenvalue problem","PEPSetProblemType",&flg));
106: if (flg) PetscCall(PEPSetProblemType(pep,PEP_HERMITIAN));
107: PetscCall(PetscOptionsBoolGroup("-pep_hyperbolic","Hyperbolic polynomial eigenvalue problem","PEPSetProblemType",&flg));
108: if (flg) PetscCall(PEPSetProblemType(pep,PEP_HYPERBOLIC));
109: PetscCall(PetscOptionsBoolGroupEnd("-pep_gyroscopic","Gyroscopic polynomial eigenvalue problem","PEPSetProblemType",&flg));
110: if (flg) PetscCall(PEPSetProblemType(pep,PEP_GYROSCOPIC));
112: scale = pep->scale;
113: PetscCall(PetscOptionsEnum("-pep_scale","Scaling strategy","PEPSetScale",PEPScaleTypes,(PetscEnum)scale,(PetscEnum*)&scale,&flg1));
114: r = pep->sfactor;
115: PetscCall(PetscOptionsReal("-pep_scale_factor","Scale factor","PEPSetScale",pep->sfactor,&r,&flg2));
116: if (!flg2 && r==1.0) r = PETSC_DEFAULT;
117: j = pep->sits;
118: PetscCall(PetscOptionsInt("-pep_scale_its","Number of iterations in diagonal scaling","PEPSetScale",pep->sits,&j,&flg3));
119: t = pep->slambda;
120: PetscCall(PetscOptionsReal("-pep_scale_lambda","Estimate of eigenvalue (modulus) for diagonal scaling","PEPSetScale",pep->slambda,&t,&flg4));
121: if (flg1 || flg2 || flg3 || flg4) PetscCall(PEPSetScale(pep,scale,r,NULL,NULL,j,t));
123: PetscCall(PetscOptionsEnum("-pep_extract","Extraction method","PEPSetExtract",PEPExtractTypes,(PetscEnum)pep->extract,(PetscEnum*)&pep->extract,NULL));
125: refine = pep->refine;
126: PetscCall(PetscOptionsEnum("-pep_refine","Iterative refinement method","PEPSetRefine",PEPRefineTypes,(PetscEnum)refine,(PetscEnum*)&refine,&flg1));
127: i = pep->npart;
128: PetscCall(PetscOptionsInt("-pep_refine_partitions","Number of partitions of the communicator for iterative refinement","PEPSetRefine",pep->npart,&i,&flg2));
129: r = pep->rtol;
130: PetscCall(PetscOptionsReal("-pep_refine_tol","Tolerance for iterative refinement","PEPSetRefine",pep->rtol==(PetscReal)PETSC_DEFAULT?SLEPC_DEFAULT_TOL/1000:pep->rtol,&r,&flg3));
131: j = pep->rits;
132: PetscCall(PetscOptionsInt("-pep_refine_its","Maximum number of iterations for iterative refinement","PEPSetRefine",pep->rits,&j,&flg4));
133: scheme = pep->scheme;
134: PetscCall(PetscOptionsEnum("-pep_refine_scheme","Scheme used for linear systems within iterative refinement","PEPSetRefine",PEPRefineSchemes,(PetscEnum)scheme,(PetscEnum*)&scheme,&flg5));
135: if (flg1 || flg2 || flg3 || flg4 || flg5) PetscCall(PEPSetRefine(pep,refine,i,r,j,scheme));
137: i = pep->max_it;
138: PetscCall(PetscOptionsInt("-pep_max_it","Maximum number of iterations","PEPSetTolerances",pep->max_it,&i,&flg1));
139: r = pep->tol;
140: PetscCall(PetscOptionsReal("-pep_tol","Tolerance","PEPSetTolerances",SlepcDefaultTol(pep->tol),&r,&flg2));
141: if (flg1 || flg2) PetscCall(PEPSetTolerances(pep,r,i));
143: PetscCall(PetscOptionsBoolGroupBegin("-pep_conv_rel","Relative error convergence test","PEPSetConvergenceTest",&flg));
144: if (flg) PetscCall(PEPSetConvergenceTest(pep,PEP_CONV_REL));
145: PetscCall(PetscOptionsBoolGroup("-pep_conv_norm","Convergence test relative to the matrix norms","PEPSetConvergenceTest",&flg));
146: if (flg) PetscCall(PEPSetConvergenceTest(pep,PEP_CONV_NORM));
147: PetscCall(PetscOptionsBoolGroup("-pep_conv_abs","Absolute error convergence test","PEPSetConvergenceTest",&flg));
148: if (flg) PetscCall(PEPSetConvergenceTest(pep,PEP_CONV_ABS));
149: PetscCall(PetscOptionsBoolGroupEnd("-pep_conv_user","User-defined convergence test","PEPSetConvergenceTest",&flg));
150: if (flg) PetscCall(PEPSetConvergenceTest(pep,PEP_CONV_USER));
152: PetscCall(PetscOptionsBoolGroupBegin("-pep_stop_basic","Stop iteration if all eigenvalues converged or max_it reached","PEPSetStoppingTest",&flg));
153: if (flg) PetscCall(PEPSetStoppingTest(pep,PEP_STOP_BASIC));
154: PetscCall(PetscOptionsBoolGroupEnd("-pep_stop_user","User-defined stopping test","PEPSetStoppingTest",&flg));
155: if (flg) PetscCall(PEPSetStoppingTest(pep,PEP_STOP_USER));
157: i = pep->nev;
158: PetscCall(PetscOptionsInt("-pep_nev","Number of eigenvalues to compute","PEPSetDimensions",pep->nev,&i,&flg1));
159: j = pep->ncv;
160: PetscCall(PetscOptionsInt("-pep_ncv","Number of basis vectors","PEPSetDimensions",pep->ncv,&j,&flg2));
161: k = pep->mpd;
162: PetscCall(PetscOptionsInt("-pep_mpd","Maximum dimension of projected problem","PEPSetDimensions",pep->mpd,&k,&flg3));
163: if (flg1 || flg2 || flg3) PetscCall(PEPSetDimensions(pep,i,j,k));
165: PetscCall(PetscOptionsEnum("-pep_basis","Polynomial basis","PEPSetBasis",PEPBasisTypes,(PetscEnum)pep->basis,(PetscEnum*)&pep->basis,NULL));
167: PetscCall(PetscOptionsBoolGroupBegin("-pep_largest_magnitude","Compute largest eigenvalues in magnitude","PEPSetWhichEigenpairs",&flg));
168: if (flg) PetscCall(PEPSetWhichEigenpairs(pep,PEP_LARGEST_MAGNITUDE));
169: PetscCall(PetscOptionsBoolGroup("-pep_smallest_magnitude","Compute smallest eigenvalues in magnitude","PEPSetWhichEigenpairs",&flg));
170: if (flg) PetscCall(PEPSetWhichEigenpairs(pep,PEP_SMALLEST_MAGNITUDE));
171: PetscCall(PetscOptionsBoolGroup("-pep_largest_real","Compute eigenvalues with largest real parts","PEPSetWhichEigenpairs",&flg));
172: if (flg) PetscCall(PEPSetWhichEigenpairs(pep,PEP_LARGEST_REAL));
173: PetscCall(PetscOptionsBoolGroup("-pep_smallest_real","Compute eigenvalues with smallest real parts","PEPSetWhichEigenpairs",&flg));
174: if (flg) PetscCall(PEPSetWhichEigenpairs(pep,PEP_SMALLEST_REAL));
175: PetscCall(PetscOptionsBoolGroup("-pep_largest_imaginary","Compute eigenvalues with largest imaginary parts","PEPSetWhichEigenpairs",&flg));
176: if (flg) PetscCall(PEPSetWhichEigenpairs(pep,PEP_LARGEST_IMAGINARY));
177: PetscCall(PetscOptionsBoolGroup("-pep_smallest_imaginary","Compute eigenvalues with smallest imaginary parts","PEPSetWhichEigenpairs",&flg));
178: if (flg) PetscCall(PEPSetWhichEigenpairs(pep,PEP_SMALLEST_IMAGINARY));
179: PetscCall(PetscOptionsBoolGroup("-pep_target_magnitude","Compute eigenvalues closest to target","PEPSetWhichEigenpairs",&flg));
180: if (flg) PetscCall(PEPSetWhichEigenpairs(pep,PEP_TARGET_MAGNITUDE));
181: PetscCall(PetscOptionsBoolGroup("-pep_target_real","Compute eigenvalues with real parts closest to target","PEPSetWhichEigenpairs",&flg));
182: if (flg) PetscCall(PEPSetWhichEigenpairs(pep,PEP_TARGET_REAL));
183: PetscCall(PetscOptionsBoolGroup("-pep_target_imaginary","Compute eigenvalues with imaginary parts closest to target","PEPSetWhichEigenpairs",&flg));
184: if (flg) PetscCall(PEPSetWhichEigenpairs(pep,PEP_TARGET_IMAGINARY));
185: PetscCall(PetscOptionsBoolGroupEnd("-pep_all","Compute all eigenvalues in an interval or a region","PEPSetWhichEigenpairs",&flg));
186: if (flg) PetscCall(PEPSetWhichEigenpairs(pep,PEP_ALL));
188: PetscCall(PetscOptionsScalar("-pep_target","Value of the target","PEPSetTarget",pep->target,&s,&flg));
189: if (flg) {
190: if (pep->which!=PEP_TARGET_REAL && pep->which!=PEP_TARGET_IMAGINARY) PetscCall(PEPSetWhichEigenpairs(pep,PEP_TARGET_MAGNITUDE));
191: PetscCall(PEPSetTarget(pep,s));
192: }
194: k = 2;
195: PetscCall(PetscOptionsRealArray("-pep_interval","Computational interval (two real values separated with a comma without spaces)","PEPSetInterval",array,&k,&flg));
196: if (flg) {
197: PetscCheck(k>1,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_SIZ,"Must pass two values in -pep_interval (comma-separated without spaces)");
198: PetscCall(PEPSetWhichEigenpairs(pep,PEP_ALL));
199: PetscCall(PEPSetInterval(pep,array[0],array[1]));
200: }
202: /* -----------------------------------------------------------------------*/
203: /*
204: Cancels all monitors hardwired into code before call to PEPSetFromOptions()
205: */
206: PetscCall(PetscOptionsBool("-pep_monitor_cancel","Remove any hardwired monitor routines","PEPMonitorCancel",PETSC_FALSE,&flg,&set));
207: if (set && flg) PetscCall(PEPMonitorCancel(pep));
208: PetscCall(PEPMonitorSetFromOptions(pep,"-pep_monitor","first_approximation",NULL,PETSC_FALSE));
209: PetscCall(PEPMonitorSetFromOptions(pep,"-pep_monitor_all","all_approximations",NULL,PETSC_TRUE));
210: PetscCall(PEPMonitorSetFromOptions(pep,"-pep_monitor_conv","convergence_history",NULL,PETSC_FALSE));
212: /* -----------------------------------------------------------------------*/
213: PetscCall(PetscOptionsName("-pep_view","Print detailed information on solver used","PEPView",&set));
214: PetscCall(PetscOptionsName("-pep_view_vectors","View computed eigenvectors","PEPVectorsView",&set));
215: PetscCall(PetscOptionsName("-pep_view_values","View computed eigenvalues","PEPValuesView",&set));
216: PetscCall(PetscOptionsName("-pep_converged_reason","Print reason for convergence, and number of iterations","PEPConvergedReasonView",&set));
217: PetscCall(PetscOptionsName("-pep_error_absolute","Print absolute errors of each eigenpair","PEPErrorView",&set));
218: PetscCall(PetscOptionsName("-pep_error_relative","Print relative errors of each eigenpair","PEPErrorView",&set));
219: PetscCall(PetscOptionsName("-pep_error_backward","Print backward errors of each eigenpair","PEPErrorView",&set));
221: PetscTryTypeMethod(pep,setfromoptions,PetscOptionsObject);
222: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)pep,PetscOptionsObject));
223: PetscOptionsEnd();
225: if (!pep->V) PetscCall(PEPGetBV(pep,&pep->V));
226: PetscCall(BVSetFromOptions(pep->V));
227: if (!pep->rg) PetscCall(PEPGetRG(pep,&pep->rg));
228: PetscCall(RGSetFromOptions(pep->rg));
229: if (!pep->ds) PetscCall(PEPGetDS(pep,&pep->ds));
230: PetscCall(PEPSetDSType(pep));
231: PetscCall(DSSetFromOptions(pep->ds));
232: if (!pep->st) PetscCall(PEPGetST(pep,&pep->st));
233: PetscCall(PEPSetDefaultST(pep));
234: PetscCall(STSetFromOptions(pep->st));
235: if (!pep->refineksp) PetscCall(PEPRefineGetKSP(pep,&pep->refineksp));
236: PetscCall(KSPSetFromOptions(pep->refineksp));
237: PetscFunctionReturn(PETSC_SUCCESS);
238: }
240: /*@C
241: PEPGetTolerances - Gets the tolerance and maximum iteration count used
242: by the PEP convergence tests.
244: Not Collective
246: Input Parameter:
247: . pep - the polynomial eigensolver context
249: Output Parameters:
250: + tol - the convergence tolerance
251: - maxits - maximum number of iterations
253: Notes:
254: The user can specify NULL for any parameter that is not needed.
256: Level: intermediate
258: .seealso: PEPSetTolerances()
259: @*/
260: PetscErrorCode PEPGetTolerances(PEP pep,PetscReal *tol,PetscInt *maxits)
261: {
262: PetscFunctionBegin;
264: if (tol) *tol = pep->tol;
265: if (maxits) *maxits = pep->max_it;
266: PetscFunctionReturn(PETSC_SUCCESS);
267: }
269: /*@
270: PEPSetTolerances - Sets the tolerance and maximum iteration count used
271: by the PEP convergence tests.
273: Logically Collective
275: Input Parameters:
276: + pep - the polynomial eigensolver context
277: . tol - the convergence tolerance
278: - maxits - maximum number of iterations to use
280: Options Database Keys:
281: + -pep_tol <tol> - Sets the convergence tolerance
282: - -pep_max_it <maxits> - Sets the maximum number of iterations allowed
284: Notes:
285: Use PETSC_DEFAULT for either argument to assign a reasonably good value.
287: Level: intermediate
289: .seealso: PEPGetTolerances()
290: @*/
291: PetscErrorCode PEPSetTolerances(PEP pep,PetscReal tol,PetscInt maxits)
292: {
293: PetscFunctionBegin;
297: if (tol == (PetscReal)PETSC_DEFAULT) {
298: pep->tol = PETSC_DEFAULT;
299: pep->state = PEP_STATE_INITIAL;
300: } else {
301: PetscCheck(tol>0.0,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
302: pep->tol = tol;
303: }
304: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
305: pep->max_it = PETSC_DEFAULT;
306: pep->state = PEP_STATE_INITIAL;
307: } else {
308: PetscCheck(maxits>0,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
309: pep->max_it = maxits;
310: }
311: PetscFunctionReturn(PETSC_SUCCESS);
312: }
314: /*@C
315: PEPGetDimensions - Gets the number of eigenvalues to compute
316: and the dimension of the subspace.
318: Not Collective
320: Input Parameter:
321: . pep - the polynomial eigensolver context
323: Output Parameters:
324: + nev - number of eigenvalues to compute
325: . ncv - the maximum dimension of the subspace to be used by the solver
326: - mpd - the maximum dimension allowed for the projected problem
328: Notes:
329: The user can specify NULL for any parameter that is not needed.
331: Level: intermediate
333: .seealso: PEPSetDimensions()
334: @*/
335: PetscErrorCode PEPGetDimensions(PEP pep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
336: {
337: PetscFunctionBegin;
339: if (nev) *nev = pep->nev;
340: if (ncv) *ncv = pep->ncv;
341: if (mpd) *mpd = pep->mpd;
342: PetscFunctionReturn(PETSC_SUCCESS);
343: }
345: /*@
346: PEPSetDimensions - Sets the number of eigenvalues to compute
347: and the dimension of the subspace.
349: Logically Collective
351: Input Parameters:
352: + pep - the polynomial eigensolver context
353: . nev - number of eigenvalues to compute
354: . ncv - the maximum dimension of the subspace to be used by the solver
355: - mpd - the maximum dimension allowed for the projected problem
357: Options Database Keys:
358: + -pep_nev <nev> - Sets the number of eigenvalues
359: . -pep_ncv <ncv> - Sets the dimension of the subspace
360: - -pep_mpd <mpd> - Sets the maximum projected dimension
362: Notes:
363: Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
364: dependent on the solution method.
366: The parameters ncv and mpd are intimately related, so that the user is advised
367: to set one of them at most. Normal usage is that
368: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
369: (b) in cases where nev is large, the user sets mpd.
371: The value of ncv should always be between nev and (nev+mpd), typically
372: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
373: a smaller value should be used.
375: When computing all eigenvalues in an interval, see PEPSetInterval(), these
376: parameters lose relevance, and tuning must be done with PEPSTOARSetDimensions().
378: Level: intermediate
380: .seealso: PEPGetDimensions(), PEPSetInterval(), PEPSTOARSetDimensions()
381: @*/
382: PetscErrorCode PEPSetDimensions(PEP pep,PetscInt nev,PetscInt ncv,PetscInt mpd)
383: {
384: PetscFunctionBegin;
389: PetscCheck(nev>0,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
390: pep->nev = nev;
391: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
392: pep->ncv = PETSC_DEFAULT;
393: } else {
394: PetscCheck(ncv>0,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
395: pep->ncv = ncv;
396: }
397: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
398: pep->mpd = PETSC_DEFAULT;
399: } else {
400: PetscCheck(mpd>0,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
401: pep->mpd = mpd;
402: }
403: pep->state = PEP_STATE_INITIAL;
404: PetscFunctionReturn(PETSC_SUCCESS);
405: }
407: /*@
408: PEPSetWhichEigenpairs - Specifies which portion of the spectrum is
409: to be sought.
411: Logically Collective
413: Input Parameters:
414: + pep - eigensolver context obtained from PEPCreate()
415: - which - the portion of the spectrum to be sought
417: Options Database Keys:
418: + -pep_largest_magnitude - Sets largest eigenvalues in magnitude
419: . -pep_smallest_magnitude - Sets smallest eigenvalues in magnitude
420: . -pep_largest_real - Sets largest real parts
421: . -pep_smallest_real - Sets smallest real parts
422: . -pep_largest_imaginary - Sets largest imaginary parts
423: . -pep_smallest_imaginary - Sets smallest imaginary parts
424: . -pep_target_magnitude - Sets eigenvalues closest to target
425: . -pep_target_real - Sets real parts closest to target
426: . -pep_target_imaginary - Sets imaginary parts closest to target
427: - -pep_all - Sets all eigenvalues in an interval or region
429: Notes:
430: The parameter 'which' can have one of these values
432: + PEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
433: . PEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
434: . PEP_LARGEST_REAL - largest real parts
435: . PEP_SMALLEST_REAL - smallest real parts
436: . PEP_LARGEST_IMAGINARY - largest imaginary parts
437: . PEP_SMALLEST_IMAGINARY - smallest imaginary parts
438: . PEP_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
439: . PEP_TARGET_REAL - eigenvalues with real part closest to target
440: . PEP_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
441: . PEP_ALL - all eigenvalues contained in a given interval or region
442: - PEP_WHICH_USER - user defined ordering set with PEPSetEigenvalueComparison()
444: Not all eigensolvers implemented in PEP account for all the possible values
445: stated above. If SLEPc is compiled for real numbers PEP_LARGEST_IMAGINARY
446: and PEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
447: for eigenvalue selection.
449: The target is a scalar value provided with PEPSetTarget().
451: The criterion PEP_TARGET_IMAGINARY is available only in case PETSc and
452: SLEPc have been built with complex scalars.
454: PEP_ALL is intended for use in combination with an interval (see
455: PEPSetInterval()), when all eigenvalues within the interval are requested,
456: and also for computing all eigenvalues in a region with the CISS solver.
457: In both cases, the number of eigenvalues is unknown, so the nev parameter
458: has a different sense, see PEPSetDimensions().
460: Level: intermediate
462: .seealso: PEPGetWhichEigenpairs(), PEPSetTarget(), PEPSetInterval(),
463: PEPSetDimensions(), PEPSetEigenvalueComparison(), PEPWhich
464: @*/
465: PetscErrorCode PEPSetWhichEigenpairs(PEP pep,PEPWhich which)
466: {
467: PetscFunctionBegin;
470: switch (which) {
471: case PEP_LARGEST_MAGNITUDE:
472: case PEP_SMALLEST_MAGNITUDE:
473: case PEP_LARGEST_REAL:
474: case PEP_SMALLEST_REAL:
475: case PEP_LARGEST_IMAGINARY:
476: case PEP_SMALLEST_IMAGINARY:
477: case PEP_TARGET_MAGNITUDE:
478: case PEP_TARGET_REAL:
479: #if defined(PETSC_USE_COMPLEX)
480: case PEP_TARGET_IMAGINARY:
481: #endif
482: case PEP_ALL:
483: case PEP_WHICH_USER:
484: if (pep->which != which) {
485: pep->state = PEP_STATE_INITIAL;
486: pep->which = which;
487: }
488: break;
489: #if !defined(PETSC_USE_COMPLEX)
490: case PEP_TARGET_IMAGINARY:
491: SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"PEP_TARGET_IMAGINARY can be used only with complex scalars");
492: #endif
493: default:
494: SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
495: }
496: PetscFunctionReturn(PETSC_SUCCESS);
497: }
499: /*@
500: PEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
501: sought.
503: Not Collective
505: Input Parameter:
506: . pep - eigensolver context obtained from PEPCreate()
508: Output Parameter:
509: . which - the portion of the spectrum to be sought
511: Notes:
512: See PEPSetWhichEigenpairs() for possible values of 'which'.
514: Level: intermediate
516: .seealso: PEPSetWhichEigenpairs(), PEPWhich
517: @*/
518: PetscErrorCode PEPGetWhichEigenpairs(PEP pep,PEPWhich *which)
519: {
520: PetscFunctionBegin;
522: PetscAssertPointer(which,2);
523: *which = pep->which;
524: PetscFunctionReturn(PETSC_SUCCESS);
525: }
527: /*@C
528: PEPSetEigenvalueComparison - Specifies the eigenvalue comparison function
529: when PEPSetWhichEigenpairs() is set to PEP_WHICH_USER.
531: Logically Collective
533: Input Parameters:
534: + pep - eigensolver context obtained from PEPCreate()
535: . comp - a pointer to the comparison function
536: - ctx - a context pointer (the last parameter to the comparison function)
538: Calling sequence of comp:
539: $ PetscErrorCode comp(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
540: + ar - real part of the 1st eigenvalue
541: . ai - imaginary part of the 1st eigenvalue
542: . br - real part of the 2nd eigenvalue
543: . bi - imaginary part of the 2nd eigenvalue
544: . res - result of comparison
545: - ctx - optional context, as set by PEPSetEigenvalueComparison()
547: Note:
548: The returning parameter 'res' can be
549: + negative - if the 1st eigenvalue is preferred to the 2st one
550: . zero - if both eigenvalues are equally preferred
551: - positive - if the 2st eigenvalue is preferred to the 1st one
553: Level: advanced
555: .seealso: PEPSetWhichEigenpairs(), PEPWhich
556: @*/
557: PetscErrorCode PEPSetEigenvalueComparison(PEP pep,PetscErrorCode (*comp)(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx),void* ctx)
558: {
559: PetscFunctionBegin;
561: pep->sc->comparison = comp;
562: pep->sc->comparisonctx = ctx;
563: pep->which = PEP_WHICH_USER;
564: PetscFunctionReturn(PETSC_SUCCESS);
565: }
567: /*@
568: PEPSetProblemType - Specifies the type of the polynomial eigenvalue problem.
570: Logically Collective
572: Input Parameters:
573: + pep - the polynomial eigensolver context
574: - type - a known type of polynomial eigenvalue problem
576: Options Database Keys:
577: + -pep_general - general problem with no particular structure
578: . -pep_hermitian - problem whose coefficient matrices are Hermitian
579: . -pep_hyperbolic - Hermitian problem that satisfies the definition of hyperbolic
580: - -pep_gyroscopic - problem with Hamiltonian structure
582: Notes:
583: Allowed values for the problem type are general (PEP_GENERAL), Hermitian
584: (PEP_HERMITIAN), hyperbolic (PEP_HYPERBOLIC), and gyroscopic (PEP_GYROSCOPIC).
586: This function is used to instruct SLEPc to exploit certain structure in
587: the polynomial eigenproblem. By default, no particular structure is assumed.
589: If the problem matrices are Hermitian (symmetric in the real case) or
590: Hermitian/skew-Hermitian then the solver can exploit this fact to perform
591: less operations or provide better stability. Hyperbolic problems are a
592: particular case of Hermitian problems, some solvers may treat them simply as
593: Hermitian.
595: Level: intermediate
597: .seealso: PEPSetOperators(), PEPSetType(), PEPGetProblemType(), PEPProblemType
598: @*/
599: PetscErrorCode PEPSetProblemType(PEP pep,PEPProblemType type)
600: {
601: PetscFunctionBegin;
604: PetscCheck(type==PEP_GENERAL || type==PEP_HERMITIAN || type==PEP_HYPERBOLIC || type==PEP_GYROSCOPIC,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
605: if (type != pep->problem_type) {
606: pep->problem_type = type;
607: pep->state = PEP_STATE_INITIAL;
608: }
609: PetscFunctionReturn(PETSC_SUCCESS);
610: }
612: /*@
613: PEPGetProblemType - Gets the problem type from the PEP object.
615: Not Collective
617: Input Parameter:
618: . pep - the polynomial eigensolver context
620: Output Parameter:
621: . type - the problem type
623: Level: intermediate
625: .seealso: PEPSetProblemType(), PEPProblemType
626: @*/
627: PetscErrorCode PEPGetProblemType(PEP pep,PEPProblemType *type)
628: {
629: PetscFunctionBegin;
631: PetscAssertPointer(type,2);
632: *type = pep->problem_type;
633: PetscFunctionReturn(PETSC_SUCCESS);
634: }
636: /*@
637: PEPSetBasis - Specifies the type of polynomial basis used to describe the
638: polynomial eigenvalue problem.
640: Logically Collective
642: Input Parameters:
643: + pep - the polynomial eigensolver context
644: - basis - the type of polynomial basis
646: Options Database Key:
647: . -pep_basis <basis> - Select the basis type
649: Notes:
650: By default, the coefficient matrices passed via PEPSetOperators() are
651: expressed in the monomial basis, i.e.
652: P(lambda) = A_0 + lambda*A_1 + lambda^2*A_2 + ... + lambda^d*A_d.
653: Other polynomial bases may have better numerical behaviour, but the user
654: must then pass the coefficient matrices accordingly.
656: Level: intermediate
658: .seealso: PEPSetOperators(), PEPGetBasis(), PEPBasis
659: @*/
660: PetscErrorCode PEPSetBasis(PEP pep,PEPBasis basis)
661: {
662: PetscFunctionBegin;
665: pep->basis = basis;
666: PetscFunctionReturn(PETSC_SUCCESS);
667: }
669: /*@
670: PEPGetBasis - Gets the type of polynomial basis from the PEP object.
672: Not Collective
674: Input Parameter:
675: . pep - the polynomial eigensolver context
677: Output Parameter:
678: . basis - the polynomial basis
680: Level: intermediate
682: .seealso: PEPSetBasis(), PEPBasis
683: @*/
684: PetscErrorCode PEPGetBasis(PEP pep,PEPBasis *basis)
685: {
686: PetscFunctionBegin;
688: PetscAssertPointer(basis,2);
689: *basis = pep->basis;
690: PetscFunctionReturn(PETSC_SUCCESS);
691: }
693: /*@
694: PEPSetTrackAll - Specifies if the solver must compute the residual of all
695: approximate eigenpairs or not.
697: Logically Collective
699: Input Parameters:
700: + pep - the eigensolver context
701: - trackall - whether compute all residuals or not
703: Notes:
704: If the user sets trackall=PETSC_TRUE then the solver explicitly computes
705: the residual for each eigenpair approximation. Computing the residual is
706: usually an expensive operation and solvers commonly compute the associated
707: residual to the first unconverged eigenpair.
709: The option '-pep_monitor_all' automatically activates this option.
711: Level: developer
713: .seealso: PEPGetTrackAll()
714: @*/
715: PetscErrorCode PEPSetTrackAll(PEP pep,PetscBool trackall)
716: {
717: PetscFunctionBegin;
720: pep->trackall = trackall;
721: PetscFunctionReturn(PETSC_SUCCESS);
722: }
724: /*@
725: PEPGetTrackAll - Returns the flag indicating whether all residual norms must
726: be computed or not.
728: Not Collective
730: Input Parameter:
731: . pep - the eigensolver context
733: Output Parameter:
734: . trackall - the returned flag
736: Level: developer
738: .seealso: PEPSetTrackAll()
739: @*/
740: PetscErrorCode PEPGetTrackAll(PEP pep,PetscBool *trackall)
741: {
742: PetscFunctionBegin;
744: PetscAssertPointer(trackall,2);
745: *trackall = pep->trackall;
746: PetscFunctionReturn(PETSC_SUCCESS);
747: }
749: /*@C
750: PEPSetConvergenceTestFunction - Sets a function to compute the error estimate
751: used in the convergence test.
753: Logically Collective
755: Input Parameters:
756: + pep - eigensolver context obtained from PEPCreate()
757: . conv - a pointer to the convergence test function
758: . ctx - context for private data for the convergence routine (may be null)
759: - destroy - a routine for destroying the context (may be null)
761: Calling sequence of conv:
762: $ PetscErrorCode conv(PEP pep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
763: + pep - eigensolver context obtained from PEPCreate()
764: . eigr - real part of the eigenvalue
765: . eigi - imaginary part of the eigenvalue
766: . res - residual norm associated to the eigenpair
767: . errest - (output) computed error estimate
768: - ctx - optional context, as set by PEPSetConvergenceTestFunction()
770: Note:
771: If the error estimate returned by the convergence test function is less than
772: the tolerance, then the eigenvalue is accepted as converged.
774: Level: advanced
776: .seealso: PEPSetConvergenceTest(), PEPSetTolerances()
777: @*/
778: PetscErrorCode PEPSetConvergenceTestFunction(PEP pep,PetscErrorCode (*conv)(PEP pep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx),void* ctx,PetscErrorCode (*destroy)(void*))
779: {
780: PetscFunctionBegin;
782: if (pep->convergeddestroy) PetscCall((*pep->convergeddestroy)(pep->convergedctx));
783: pep->convergeduser = conv;
784: pep->convergeddestroy = destroy;
785: pep->convergedctx = ctx;
786: if (conv == PEPConvergedRelative) pep->conv = PEP_CONV_REL;
787: else if (conv == PEPConvergedNorm) pep->conv = PEP_CONV_NORM;
788: else if (conv == PEPConvergedAbsolute) pep->conv = PEP_CONV_ABS;
789: else {
790: pep->conv = PEP_CONV_USER;
791: pep->converged = pep->convergeduser;
792: }
793: PetscFunctionReturn(PETSC_SUCCESS);
794: }
796: /*@
797: PEPSetConvergenceTest - Specifies how to compute the error estimate
798: used in the convergence test.
800: Logically Collective
802: Input Parameters:
803: + pep - eigensolver context obtained from PEPCreate()
804: - conv - the type of convergence test
806: Options Database Keys:
807: + -pep_conv_abs - Sets the absolute convergence test
808: . -pep_conv_rel - Sets the convergence test relative to the eigenvalue
809: . -pep_conv_norm - Sets the convergence test relative to the matrix norms
810: - -pep_conv_user - Selects the user-defined convergence test
812: Note:
813: The parameter 'conv' can have one of these values
814: + PEP_CONV_ABS - absolute error ||r||
815: . PEP_CONV_REL - error relative to the eigenvalue l, ||r||/|l|
816: . PEP_CONV_NORM - error relative matrix norms, ||r||/sum_i(l^i*||A_i||)
817: - PEP_CONV_USER - function set by PEPSetConvergenceTestFunction()
819: Level: intermediate
821: .seealso: PEPGetConvergenceTest(), PEPSetConvergenceTestFunction(), PEPSetStoppingTest(), PEPConv
822: @*/
823: PetscErrorCode PEPSetConvergenceTest(PEP pep,PEPConv conv)
824: {
825: PetscFunctionBegin;
828: switch (conv) {
829: case PEP_CONV_ABS: pep->converged = PEPConvergedAbsolute; break;
830: case PEP_CONV_REL: pep->converged = PEPConvergedRelative; break;
831: case PEP_CONV_NORM: pep->converged = PEPConvergedNorm; break;
832: case PEP_CONV_USER:
833: PetscCheck(pep->convergeduser,PetscObjectComm((PetscObject)pep),PETSC_ERR_ORDER,"Must call PEPSetConvergenceTestFunction() first");
834: pep->converged = pep->convergeduser;
835: break;
836: default:
837: SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
838: }
839: pep->conv = conv;
840: PetscFunctionReturn(PETSC_SUCCESS);
841: }
843: /*@
844: PEPGetConvergenceTest - Gets the method used to compute the error estimate
845: used in the convergence test.
847: Not Collective
849: Input Parameters:
850: . pep - eigensolver context obtained from PEPCreate()
852: Output Parameters:
853: . conv - the type of convergence test
855: Level: intermediate
857: .seealso: PEPSetConvergenceTest(), PEPConv
858: @*/
859: PetscErrorCode PEPGetConvergenceTest(PEP pep,PEPConv *conv)
860: {
861: PetscFunctionBegin;
863: PetscAssertPointer(conv,2);
864: *conv = pep->conv;
865: PetscFunctionReturn(PETSC_SUCCESS);
866: }
868: /*@C
869: PEPSetStoppingTestFunction - Sets a function to decide when to stop the outer
870: iteration of the eigensolver.
872: Logically Collective
874: Input Parameters:
875: + pep - eigensolver context obtained from PEPCreate()
876: . stop - pointer to the stopping test function
877: . ctx - context for private data for the stopping routine (may be null)
878: - destroy - a routine for destroying the context (may be null)
880: Calling sequence of stop:
881: $ PetscErrorCode stop(PEP pep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,PEPConvergedReason *reason,void *ctx)
882: + pep - eigensolver context obtained from PEPCreate()
883: . its - current number of iterations
884: . max_it - maximum number of iterations
885: . nconv - number of currently converged eigenpairs
886: . nev - number of requested eigenpairs
887: . reason - (output) result of the stopping test
888: - ctx - optional context, as set by PEPSetStoppingTestFunction()
890: Note:
891: Normal usage is to first call the default routine PEPStoppingBasic() and then
892: set reason to PEP_CONVERGED_USER if some user-defined conditions have been
893: met. To let the eigensolver continue iterating, the result must be left as
894: PEP_CONVERGED_ITERATING.
896: Level: advanced
898: .seealso: PEPSetStoppingTest(), PEPStoppingBasic()
899: @*/
900: PetscErrorCode PEPSetStoppingTestFunction(PEP pep,PetscErrorCode (*stop)(PEP pep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,PEPConvergedReason *reason,void *ctx),void* ctx,PetscErrorCode (*destroy)(void*))
901: {
902: PetscFunctionBegin;
904: if (pep->stoppingdestroy) PetscCall((*pep->stoppingdestroy)(pep->stoppingctx));
905: pep->stoppinguser = stop;
906: pep->stoppingdestroy = destroy;
907: pep->stoppingctx = ctx;
908: if (stop == PEPStoppingBasic) pep->stop = PEP_STOP_BASIC;
909: else {
910: pep->stop = PEP_STOP_USER;
911: pep->stopping = pep->stoppinguser;
912: }
913: PetscFunctionReturn(PETSC_SUCCESS);
914: }
916: /*@
917: PEPSetStoppingTest - Specifies how to decide the termination of the outer
918: loop of the eigensolver.
920: Logically Collective
922: Input Parameters:
923: + pep - eigensolver context obtained from PEPCreate()
924: - stop - the type of stopping test
926: Options Database Keys:
927: + -pep_stop_basic - Sets the default stopping test
928: - -pep_stop_user - Selects the user-defined stopping test
930: Note:
931: The parameter 'stop' can have one of these values
932: + PEP_STOP_BASIC - default stopping test
933: - PEP_STOP_USER - function set by PEPSetStoppingTestFunction()
935: Level: advanced
937: .seealso: PEPGetStoppingTest(), PEPSetStoppingTestFunction(), PEPSetConvergenceTest(), PEPStop
938: @*/
939: PetscErrorCode PEPSetStoppingTest(PEP pep,PEPStop stop)
940: {
941: PetscFunctionBegin;
944: switch (stop) {
945: case PEP_STOP_BASIC: pep->stopping = PEPStoppingBasic; break;
946: case PEP_STOP_USER:
947: PetscCheck(pep->stoppinguser,PetscObjectComm((PetscObject)pep),PETSC_ERR_ORDER,"Must call PEPSetStoppingTestFunction() first");
948: pep->stopping = pep->stoppinguser;
949: break;
950: default:
951: SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
952: }
953: pep->stop = stop;
954: PetscFunctionReturn(PETSC_SUCCESS);
955: }
957: /*@
958: PEPGetStoppingTest - Gets the method used to decide the termination of the outer
959: loop of the eigensolver.
961: Not Collective
963: Input Parameters:
964: . pep - eigensolver context obtained from PEPCreate()
966: Output Parameters:
967: . stop - the type of stopping test
969: Level: advanced
971: .seealso: PEPSetStoppingTest(), PEPStop
972: @*/
973: PetscErrorCode PEPGetStoppingTest(PEP pep,PEPStop *stop)
974: {
975: PetscFunctionBegin;
977: PetscAssertPointer(stop,2);
978: *stop = pep->stop;
979: PetscFunctionReturn(PETSC_SUCCESS);
980: }
982: /*@
983: PEPSetScale - Specifies the scaling strategy to be used.
985: Collective
987: Input Parameters:
988: + pep - the eigensolver context
989: . scale - scaling strategy
990: . alpha - the scaling factor used in the scalar strategy
991: . Dl - the left diagonal matrix of the diagonal scaling algorithm
992: . Dr - the right diagonal matrix of the diagonal scaling algorithm
993: . its - number of iterations of the diagonal scaling algorithm
994: - lambda - approximation to wanted eigenvalues (modulus)
996: Options Database Keys:
997: + -pep_scale <type> - scaling type, one of <none,scalar,diagonal,both>
998: . -pep_scale_factor <alpha> - the scaling factor
999: . -pep_scale_its <its> - number of iterations
1000: - -pep_scale_lambda <lambda> - approximation to eigenvalues
1002: Notes:
1003: There are two non-exclusive scaling strategies, scalar and diagonal.
1005: In the scalar strategy, scaling is applied to the eigenvalue, that is,
1006: mu = lambda/alpha is the new eigenvalue and all matrices are scaled
1007: accordingly. After solving the scaled problem, the original lambda is
1008: recovered. Parameter 'alpha' must be positive. Use PETSC_DEFAULT to let
1009: the solver compute a reasonable scaling factor.
1011: In the diagonal strategy, the solver works implicitly with matrix Dl*A*Dr,
1012: where Dl and Dr are appropriate diagonal matrices. This improves the accuracy
1013: of the computed results in some cases. The user may provide the Dr and Dl
1014: matrices represented as Vec objects storing diagonal elements. If not
1015: provided, these matrices are computed internally. This option requires
1016: that the polynomial coefficient matrices are of MATAIJ type.
1017: The parameter 'its' is the number of iterations performed by the method.
1018: Parameter 'lambda' must be positive. Use PETSC_DEFAULT or set lambda = 1.0 if
1019: no information about eigenvalues is available.
1021: Level: intermediate
1023: .seealso: PEPGetScale()
1024: @*/
1025: PetscErrorCode PEPSetScale(PEP pep,PEPScale scale,PetscReal alpha,Vec Dl,Vec Dr,PetscInt its,PetscReal lambda)
1026: {
1027: PetscFunctionBegin;
1030: pep->scale = scale;
1031: if (scale==PEP_SCALE_SCALAR || scale==PEP_SCALE_BOTH) {
1033: if (alpha == (PetscReal)PETSC_DEFAULT || alpha == (PetscReal)PETSC_DECIDE) {
1034: pep->sfactor = 0.0;
1035: pep->sfactor_set = PETSC_FALSE;
1036: } else {
1037: PetscCheck(alpha>0.0,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of alpha. Must be > 0");
1038: pep->sfactor = alpha;
1039: pep->sfactor_set = PETSC_TRUE;
1040: }
1041: }
1042: if (scale==PEP_SCALE_DIAGONAL || scale==PEP_SCALE_BOTH) {
1043: if (Dl) {
1045: PetscCheckSameComm(pep,1,Dl,4);
1046: PetscCall(PetscObjectReference((PetscObject)Dl));
1047: PetscCall(VecDestroy(&pep->Dl));
1048: pep->Dl = Dl;
1049: }
1050: if (Dr) {
1052: PetscCheckSameComm(pep,1,Dr,5);
1053: PetscCall(PetscObjectReference((PetscObject)Dr));
1054: PetscCall(VecDestroy(&pep->Dr));
1055: pep->Dr = Dr;
1056: }
1059: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) pep->sits = 5;
1060: else pep->sits = its;
1061: if (lambda == (PetscReal)PETSC_DECIDE || lambda == (PetscReal)PETSC_DEFAULT) pep->slambda = 1.0;
1062: else {
1063: PetscCheck(lambda>0.0,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of lambda. Must be > 0");
1064: pep->slambda = lambda;
1065: }
1066: }
1067: pep->state = PEP_STATE_INITIAL;
1068: PetscFunctionReturn(PETSC_SUCCESS);
1069: }
1071: /*@C
1072: PEPGetScale - Gets the scaling strategy used by the PEP object, and the
1073: associated parameters.
1075: Not Collective
1077: Input Parameter:
1078: . pep - the eigensolver context
1080: Output Parameters:
1081: + scale - scaling strategy
1082: . alpha - the scaling factor used in the scalar strategy
1083: . Dl - the left diagonal matrix of the diagonal scaling algorithm
1084: . Dr - the right diagonal matrix of the diagonal scaling algorithm
1085: . its - number of iterations of the diagonal scaling algorithm
1086: - lambda - approximation to wanted eigenvalues (modulus)
1088: Level: intermediate
1090: Note:
1091: The user can specify NULL for any parameter that is not needed.
1093: If Dl or Dr were not set by the user, then the ones computed internally are
1094: returned (or a null pointer if called before PEPSetUp).
1096: .seealso: PEPSetScale(), PEPSetUp()
1097: @*/
1098: PetscErrorCode PEPGetScale(PEP pep,PEPScale *scale,PetscReal *alpha,Vec *Dl,Vec *Dr,PetscInt *its,PetscReal *lambda)
1099: {
1100: PetscFunctionBegin;
1102: if (scale) *scale = pep->scale;
1103: if (alpha) *alpha = pep->sfactor;
1104: if (Dl) *Dl = pep->Dl;
1105: if (Dr) *Dr = pep->Dr;
1106: if (its) *its = pep->sits;
1107: if (lambda) *lambda = pep->slambda;
1108: PetscFunctionReturn(PETSC_SUCCESS);
1109: }
1111: /*@
1112: PEPSetExtract - Specifies the extraction strategy to be used.
1114: Logically Collective
1116: Input Parameters:
1117: + pep - the eigensolver context
1118: - extract - extraction strategy
1120: Options Database Keys:
1121: . -pep_extract <type> - extraction type, one of <none,norm,residual,structured>
1123: Level: intermediate
1125: .seealso: PEPGetExtract()
1126: @*/
1127: PetscErrorCode PEPSetExtract(PEP pep,PEPExtract extract)
1128: {
1129: PetscFunctionBegin;
1132: pep->extract = extract;
1133: PetscFunctionReturn(PETSC_SUCCESS);
1134: }
1136: /*@
1137: PEPGetExtract - Gets the extraction strategy used by the PEP object.
1139: Not Collective
1141: Input Parameter:
1142: . pep - the eigensolver context
1144: Output Parameter:
1145: . extract - extraction strategy
1147: Level: intermediate
1149: .seealso: PEPSetExtract(), PEPExtract
1150: @*/
1151: PetscErrorCode PEPGetExtract(PEP pep,PEPExtract *extract)
1152: {
1153: PetscFunctionBegin;
1155: PetscAssertPointer(extract,2);
1156: *extract = pep->extract;
1157: PetscFunctionReturn(PETSC_SUCCESS);
1158: }
1160: /*@
1161: PEPSetRefine - Specifies the refinement type (and options) to be used
1162: after the solve.
1164: Logically Collective
1166: Input Parameters:
1167: + pep - the polynomial eigensolver context
1168: . refine - refinement type
1169: . npart - number of partitions of the communicator
1170: . tol - the convergence tolerance
1171: . its - maximum number of refinement iterations
1172: - scheme - which scheme to be used for solving the involved linear systems
1174: Options Database Keys:
1175: + -pep_refine <type> - refinement type, one of <none,simple,multiple>
1176: . -pep_refine_partitions <n> - the number of partitions
1177: . -pep_refine_tol <tol> - the tolerance
1178: . -pep_refine_its <its> - number of iterations
1179: - -pep_refine_scheme - to set the scheme for the linear solves
1181: Notes:
1182: By default, iterative refinement is disabled, since it may be very
1183: costly. There are two possible refinement strategies, simple and multiple.
1184: The simple approach performs iterative refinement on each of the
1185: converged eigenpairs individually, whereas the multiple strategy works
1186: with the invariant pair as a whole, refining all eigenpairs simultaneously.
1187: The latter may be required for the case of multiple eigenvalues.
1189: In some cases, especially when using direct solvers within the
1190: iterative refinement method, it may be helpful for improved scalability
1191: to split the communicator in several partitions. The npart parameter
1192: indicates how many partitions to use (defaults to 1).
1194: The tol and its parameters specify the stopping criterion. In the simple
1195: method, refinement continues until the residual of each eigenpair is
1196: below the tolerance (tol defaults to the PEP tol, but may be set to a
1197: different value). In contrast, the multiple method simply performs its
1198: refinement iterations (just one by default).
1200: The scheme argument is used to change the way in which linear systems are
1201: solved. Possible choices are explicit, mixed block elimination (MBE),
1202: and Schur complement.
1204: Level: intermediate
1206: .seealso: PEPGetRefine()
1207: @*/
1208: PetscErrorCode PEPSetRefine(PEP pep,PEPRefine refine,PetscInt npart,PetscReal tol,PetscInt its,PEPRefineScheme scheme)
1209: {
1210: PetscMPIInt size;
1212: PetscFunctionBegin;
1219: pep->refine = refine;
1220: if (refine) { /* process parameters only if not REFINE_NONE */
1221: if (npart!=pep->npart) {
1222: PetscCall(PetscSubcommDestroy(&pep->refinesubc));
1223: PetscCall(KSPDestroy(&pep->refineksp));
1224: }
1225: if (npart == PETSC_DEFAULT || npart == PETSC_DECIDE) {
1226: pep->npart = 1;
1227: } else {
1228: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)pep),&size));
1229: PetscCheck(npart>0 && npart<=size,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of npart");
1230: pep->npart = npart;
1231: }
1232: if (tol == (PetscReal)PETSC_DEFAULT || tol == (PetscReal)PETSC_DECIDE) {
1233: pep->rtol = PETSC_DEFAULT;
1234: } else {
1235: PetscCheck(tol>0.0,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
1236: pep->rtol = tol;
1237: }
1238: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) {
1239: pep->rits = PETSC_DEFAULT;
1240: } else {
1241: PetscCheck(its>=0,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be >= 0");
1242: pep->rits = its;
1243: }
1244: pep->scheme = scheme;
1245: }
1246: pep->state = PEP_STATE_INITIAL;
1247: PetscFunctionReturn(PETSC_SUCCESS);
1248: }
1250: /*@C
1251: PEPGetRefine - Gets the refinement strategy used by the PEP object, and the
1252: associated parameters.
1254: Not Collective
1256: Input Parameter:
1257: . pep - the polynomial eigensolver context
1259: Output Parameters:
1260: + refine - refinement type
1261: . npart - number of partitions of the communicator
1262: . tol - the convergence tolerance
1263: . its - maximum number of refinement iterations
1264: - scheme - the scheme used for solving linear systems
1266: Level: intermediate
1268: Note:
1269: The user can specify NULL for any parameter that is not needed.
1271: .seealso: PEPSetRefine()
1272: @*/
1273: PetscErrorCode PEPGetRefine(PEP pep,PEPRefine *refine,PetscInt *npart,PetscReal *tol,PetscInt *its,PEPRefineScheme *scheme)
1274: {
1275: PetscFunctionBegin;
1277: if (refine) *refine = pep->refine;
1278: if (npart) *npart = pep->npart;
1279: if (tol) *tol = pep->rtol;
1280: if (its) *its = pep->rits;
1281: if (scheme) *scheme = pep->scheme;
1282: PetscFunctionReturn(PETSC_SUCCESS);
1283: }
1285: /*@C
1286: PEPSetOptionsPrefix - Sets the prefix used for searching for all
1287: PEP options in the database.
1289: Logically Collective
1291: Input Parameters:
1292: + pep - the polynomial eigensolver context
1293: - prefix - the prefix string to prepend to all PEP option requests
1295: Notes:
1296: A hyphen (-) must NOT be given at the beginning of the prefix name.
1297: The first character of all runtime options is AUTOMATICALLY the
1298: hyphen.
1300: For example, to distinguish between the runtime options for two
1301: different PEP contexts, one could call
1302: .vb
1303: PEPSetOptionsPrefix(pep1,"qeig1_")
1304: PEPSetOptionsPrefix(pep2,"qeig2_")
1305: .ve
1307: Level: advanced
1309: .seealso: PEPAppendOptionsPrefix(), PEPGetOptionsPrefix()
1310: @*/
1311: PetscErrorCode PEPSetOptionsPrefix(PEP pep,const char *prefix)
1312: {
1313: PetscFunctionBegin;
1315: if (!pep->st) PetscCall(PEPGetST(pep,&pep->st));
1316: PetscCall(STSetOptionsPrefix(pep->st,prefix));
1317: if (!pep->V) PetscCall(PEPGetBV(pep,&pep->V));
1318: PetscCall(BVSetOptionsPrefix(pep->V,prefix));
1319: if (!pep->ds) PetscCall(PEPGetDS(pep,&pep->ds));
1320: PetscCall(DSSetOptionsPrefix(pep->ds,prefix));
1321: if (!pep->rg) PetscCall(PEPGetRG(pep,&pep->rg));
1322: PetscCall(RGSetOptionsPrefix(pep->rg,prefix));
1323: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)pep,prefix));
1324: PetscFunctionReturn(PETSC_SUCCESS);
1325: }
1327: /*@C
1328: PEPAppendOptionsPrefix - Appends to the prefix used for searching for all
1329: PEP options in the database.
1331: Logically Collective
1333: Input Parameters:
1334: + pep - the polynomial eigensolver context
1335: - prefix - the prefix string to prepend to all PEP option requests
1337: Notes:
1338: A hyphen (-) must NOT be given at the beginning of the prefix name.
1339: The first character of all runtime options is AUTOMATICALLY the hyphen.
1341: Level: advanced
1343: .seealso: PEPSetOptionsPrefix(), PEPGetOptionsPrefix()
1344: @*/
1345: PetscErrorCode PEPAppendOptionsPrefix(PEP pep,const char *prefix)
1346: {
1347: PetscFunctionBegin;
1349: if (!pep->st) PetscCall(PEPGetST(pep,&pep->st));
1350: PetscCall(STAppendOptionsPrefix(pep->st,prefix));
1351: if (!pep->V) PetscCall(PEPGetBV(pep,&pep->V));
1352: PetscCall(BVAppendOptionsPrefix(pep->V,prefix));
1353: if (!pep->ds) PetscCall(PEPGetDS(pep,&pep->ds));
1354: PetscCall(DSAppendOptionsPrefix(pep->ds,prefix));
1355: if (!pep->rg) PetscCall(PEPGetRG(pep,&pep->rg));
1356: PetscCall(RGAppendOptionsPrefix(pep->rg,prefix));
1357: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)pep,prefix));
1358: PetscFunctionReturn(PETSC_SUCCESS);
1359: }
1361: /*@C
1362: PEPGetOptionsPrefix - Gets the prefix used for searching for all
1363: PEP options in the database.
1365: Not Collective
1367: Input Parameters:
1368: . pep - the polynomial eigensolver context
1370: Output Parameters:
1371: . prefix - pointer to the prefix string used is returned
1373: Note:
1374: On the Fortran side, the user should pass in a string 'prefix' of
1375: sufficient length to hold the prefix.
1377: Level: advanced
1379: .seealso: PEPSetOptionsPrefix(), PEPAppendOptionsPrefix()
1380: @*/
1381: PetscErrorCode PEPGetOptionsPrefix(PEP pep,const char *prefix[])
1382: {
1383: PetscFunctionBegin;
1385: PetscAssertPointer(prefix,2);
1386: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)pep,prefix));
1387: PetscFunctionReturn(PETSC_SUCCESS);
1388: }