1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2018, 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: NEP routines related to options that can be set via the command-line
12: or procedurally
13: */
15: #include <slepc/private/nepimpl.h> /*I "slepcnep.h" I*/
16: #include <petscdraw.h>
18: /*@C
19: NEPMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
20: indicated by the user.
22: Collective on NEP 24: Input Parameters:
25: + nep - the nonlinear eigensolver context
26: . name - the monitor option name
27: . help - message indicating what monitoring is done
28: . manual - manual page for the monitor
29: . monitor - the monitor function, whose context is a PetscViewerAndFormat
30: - trackall - whether this monitor tracks all eigenvalues or not
32: Level: developer
34: .seealso: NEPMonitorSet(), NEPSetTrackAll(), NEPConvMonitorSetFromOptions()
35: @*/
36: PetscErrorCode NEPMonitorSetFromOptions(NEP nep,const char name[],const char help[],const char manual[],PetscErrorCode (*monitor)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,PetscViewerAndFormat*),PetscBool trackall) 37: {
38: PetscErrorCode ierr;
39: PetscBool flg;
40: PetscViewer viewer;
41: PetscViewerFormat format;
42: PetscViewerAndFormat *vf;
45: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,name,&viewer,&format,&flg);
46: if (flg) {
47: PetscViewerAndFormatCreate(viewer,format,&vf);
48: PetscObjectDereference((PetscObject)viewer);
49: NEPMonitorSet(nep,(PetscErrorCode (*)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))monitor,vf,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
50: if (trackall) {
51: NEPSetTrackAll(nep,PETSC_TRUE);
52: }
53: }
54: return(0);
55: }
57: /*@C
58: NEPConvMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
59: indicated by the user (for monitors that only show iteration numbers of convergence).
61: Collective on NEP 63: Input Parameters:
64: + nep - the nonlinear eigensolver context
65: . name - the monitor option name
66: . help - message indicating what monitoring is done
67: . manual - manual page for the monitor
68: - monitor - the monitor function, whose context is a SlepcConvMonitor
70: Level: developer
72: .seealso: NEPMonitorSet(), NEPMonitorSetFromOptions()
73: @*/
74: PetscErrorCode NEPConvMonitorSetFromOptions(NEP nep,const char name[],const char help[],const char manual[],PetscErrorCode (*monitor)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,SlepcConvMonitor)) 75: {
76: PetscErrorCode ierr;
77: PetscBool flg;
78: PetscViewer viewer;
79: PetscViewerFormat format;
80: SlepcConvMonitor ctx;
83: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,name,&viewer,&format,&flg);
84: if (flg) {
85: SlepcConvMonitorCreate(viewer,format,&ctx);
86: PetscObjectDereference((PetscObject)viewer);
87: NEPMonitorSet(nep,(PetscErrorCode (*)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))monitor,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
88: }
89: return(0);
90: }
92: /*@
93: NEPSetFromOptions - Sets NEP options from the options database.
94: This routine must be called before NEPSetUp() if the user is to be
95: allowed to set the solver type.
97: Collective on NEP 99: Input Parameters:
100: . nep - the nonlinear eigensolver context
102: Notes:
103: To see all options, run your program with the -help option.
105: Level: beginner
106: @*/
107: PetscErrorCode NEPSetFromOptions(NEP nep)108: {
109: PetscErrorCode ierr;
110: char type[256];
111: PetscBool set,flg,flg1,flg2,flg3,flg4,flg5;
112: PetscReal r;
113: PetscScalar s;
114: PetscInt i,j,k;
115: PetscDrawLG lg;
116: NEPRefine refine;
117: NEPRefineScheme scheme;
121: NEPRegisterAll();
122: PetscObjectOptionsBegin((PetscObject)nep);
123: PetscOptionsFList("-nep_type","Nonlinear eigensolver method","NEPSetType",NEPList,(char*)(((PetscObject)nep)->type_name?((PetscObject)nep)->type_name:NEPRII),type,256,&flg);
124: if (flg) {
125: NEPSetType(nep,type);
126: } else if (!((PetscObject)nep)->type_name) {
127: NEPSetType(nep,NEPRII);
128: }
130: PetscOptionsBoolGroupBegin("-nep_general","General nonlinear eigenvalue problem","NEPSetProblemType",&flg);
131: if (flg) { NEPSetProblemType(nep,NEP_GENERAL); }
132: PetscOptionsBoolGroupEnd("-nep_rational","Rational eigenvalue problem","NEPSetProblemType",&flg);
133: if (flg) { NEPSetProblemType(nep,NEP_RATIONAL); }
135: refine = nep->refine;
136: PetscOptionsEnum("-nep_refine","Iterative refinement method","NEPSetRefine",NEPRefineTypes,(PetscEnum)refine,(PetscEnum*)&refine,&flg1);
137: i = nep->npart;
138: PetscOptionsInt("-nep_refine_partitions","Number of partitions of the communicator for iterative refinement","NEPSetRefine",nep->npart,&i,&flg2);
139: r = nep->rtol;
140: PetscOptionsReal("-nep_refine_tol","Tolerance for iterative refinement","NEPSetRefine",nep->rtol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL/1000:nep->rtol,&r,&flg3);
141: j = nep->rits;
142: PetscOptionsInt("-nep_refine_its","Maximum number of iterations for iterative refinement","NEPSetRefine",nep->rits,&j,&flg4);
143: scheme = nep->scheme;
144: PetscOptionsEnum("-nep_refine_scheme","Scheme used for linear systems within iterative refinement","NEPSetRefine",NEPRefineSchemes,(PetscEnum)scheme,(PetscEnum*)&scheme,&flg5);
145: if (flg1 || flg2 || flg3 || flg4 || flg5) { NEPSetRefine(nep,refine,i,r,j,scheme); }
147: i = nep->max_it? nep->max_it: PETSC_DEFAULT;
148: PetscOptionsInt("-nep_max_it","Maximum number of iterations","NEPSetTolerances",nep->max_it,&i,&flg1);
149: r = nep->tol;
150: PetscOptionsReal("-nep_tol","Tolerance","NEPSetTolerances",nep->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:nep->tol,&r,&flg2);
151: if (flg1 || flg2) { NEPSetTolerances(nep,r,i); }
153: PetscOptionsBoolGroupBegin("-nep_conv_rel","Relative error convergence test","NEPSetConvergenceTest",&flg);
154: if (flg) { NEPSetConvergenceTest(nep,NEP_CONV_REL); }
155: PetscOptionsBoolGroup("-nep_conv_norm","Convergence test relative to the matrix norms","NEPSetConvergenceTest",&flg);
156: if (flg) { NEPSetConvergenceTest(nep,NEP_CONV_NORM); }
157: PetscOptionsBoolGroup("-nep_conv_abs","Absolute error convergence test","NEPSetConvergenceTest",&flg);
158: if (flg) { NEPSetConvergenceTest(nep,NEP_CONV_ABS); }
159: PetscOptionsBoolGroupEnd("-nep_conv_user","User-defined convergence test","NEPSetConvergenceTest",&flg);
160: if (flg) { NEPSetConvergenceTest(nep,NEP_CONV_USER); }
162: PetscOptionsBoolGroupBegin("-nep_stop_basic","Stop iteration if all eigenvalues converged or max_it reached","NEPSetStoppingTest",&flg);
163: if (flg) { NEPSetStoppingTest(nep,NEP_STOP_BASIC); }
164: PetscOptionsBoolGroupEnd("-nep_stop_user","User-defined stopping test","NEPSetStoppingTest",&flg);
165: if (flg) { NEPSetStoppingTest(nep,NEP_STOP_USER); }
167: i = nep->nev;
168: PetscOptionsInt("-nep_nev","Number of eigenvalues to compute","NEPSetDimensions",nep->nev,&i,&flg1);
169: j = nep->ncv? nep->ncv: PETSC_DEFAULT;
170: PetscOptionsInt("-nep_ncv","Number of basis vectors","NEPSetDimensions",nep->ncv,&j,&flg2);
171: k = nep->mpd? nep->mpd: PETSC_DEFAULT;
172: PetscOptionsInt("-nep_mpd","Maximum dimension of projected problem","NEPSetDimensions",nep->mpd,&k,&flg3);
173: if (flg1 || flg2 || flg3) {
174: NEPSetDimensions(nep,i,j,k);
175: }
177: PetscOptionsBoolGroupBegin("-nep_largest_magnitude","Compute largest eigenvalues in magnitude","NEPSetWhichEigenpairs",&flg);
178: if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_MAGNITUDE); }
179: PetscOptionsBoolGroup("-nep_smallest_magnitude","Compute smallest eigenvalues in magnitude","NEPSetWhichEigenpairs",&flg);
180: if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_MAGNITUDE); }
181: PetscOptionsBoolGroup("-nep_largest_real","Compute eigenvalues with largest real parts","NEPSetWhichEigenpairs",&flg);
182: if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_REAL); }
183: PetscOptionsBoolGroup("-nep_smallest_real","Compute eigenvalues with smallest real parts","NEPSetWhichEigenpairs",&flg);
184: if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_REAL); }
185: PetscOptionsBoolGroup("-nep_largest_imaginary","Compute eigenvalues with largest imaginary parts","NEPSetWhichEigenpairs",&flg);
186: if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_IMAGINARY); }
187: PetscOptionsBoolGroup("-nep_smallest_imaginary","Compute eigenvalues with smallest imaginary parts","NEPSetWhichEigenpairs",&flg);
188: if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_IMAGINARY); }
189: PetscOptionsBoolGroup("-nep_target_magnitude","Compute eigenvalues closest to target","NEPSetWhichEigenpairs",&flg);
190: if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE); }
191: PetscOptionsBoolGroup("-nep_target_real","Compute eigenvalues with real parts closest to target","NEPSetWhichEigenpairs",&flg);
192: if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_REAL); }
193: PetscOptionsBoolGroup("-nep_target_imaginary","Compute eigenvalues with imaginary parts closest to target","NEPSetWhichEigenpairs",&flg);
194: if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_IMAGINARY); }
195: PetscOptionsBoolGroupEnd("-nep_all","Compute all eigenvalues in a region","NEPSetWhichEigenpairs",&flg);
196: if (flg) { NEPSetWhichEigenpairs(nep,NEP_ALL); }
198: PetscOptionsScalar("-nep_target","Value of the target","NEPSetTarget",nep->target,&s,&flg);
199: if (flg) {
200: if (nep->which!=NEP_TARGET_REAL && nep->which!=NEP_TARGET_IMAGINARY) {
201: NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE);
202: }
203: NEPSetTarget(nep,s);
204: }
206: /* -----------------------------------------------------------------------*/
207: /*
208: Cancels all monitors hardwired into code before call to NEPSetFromOptions()
209: */
210: PetscOptionsBool("-nep_monitor_cancel","Remove any hardwired monitor routines","NEPMonitorCancel",PETSC_FALSE,&flg,&set);
211: if (set && flg) {
212: NEPMonitorCancel(nep);
213: }
214: /*
215: Text monitors
216: */
217: NEPMonitorSetFromOptions(nep,"-nep_monitor","Monitor first unconverged approximate eigenvalue and error estimate","NEPMonitorFirst",NEPMonitorFirst,PETSC_FALSE);
218: NEPConvMonitorSetFromOptions(nep,"-nep_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","NEPMonitorConverged",NEPMonitorConverged);
219: NEPMonitorSetFromOptions(nep,"-nep_monitor_all","Monitor approximate eigenvalues and error estimates","NEPMonitorAll",NEPMonitorAll,PETSC_TRUE);
220: /*
221: Line graph monitors
222: */
223: PetscOptionsBool("-nep_monitor_lg","Monitor first unconverged approximate error estimate graphically","NEPMonitorSet",PETSC_FALSE,&flg,&set);
224: if (set && flg) {
225: NEPMonitorLGCreate(PetscObjectComm((PetscObject)nep),NULL,"Error estimates",PETSC_DECIDE,PETSC_DECIDE,300,300,&lg);
226: NEPMonitorSet(nep,NEPMonitorLG,lg,(PetscErrorCode (*)(void**))PetscDrawLGDestroy);
227: }
228: PetscOptionsBool("-nep_monitor_lg_all","Monitor error estimates graphically","NEPMonitorSet",PETSC_FALSE,&flg,&set);
229: if (set && flg) {
230: NEPMonitorLGCreate(PetscObjectComm((PetscObject)nep),NULL,"Error estimates",PETSC_DECIDE,PETSC_DECIDE,300,300,&lg);
231: NEPMonitorSet(nep,NEPMonitorLGAll,lg,(PetscErrorCode (*)(void**))PetscDrawLGDestroy);
232: NEPSetTrackAll(nep,PETSC_TRUE);
233: }
235: /* -----------------------------------------------------------------------*/
236: PetscOptionsName("-nep_view","Print detailed information on solver used","NEPView",NULL);
237: PetscOptionsName("-nep_view_vectors","View computed eigenvectors","NEPVectorsView",NULL);
238: PetscOptionsName("-nep_view_values","View computed eigenvalues","NEPValuesView",NULL);
239: PetscOptionsName("-nep_converged_reason","Print reason for convergence, and number of iterations","NEPReasonView",NULL);
240: PetscOptionsName("-nep_error_absolute","Print absolute errors of each eigenpair","NEPErrorView",NULL);
241: PetscOptionsName("-nep_error_relative","Print relative errors of each eigenpair","NEPErrorView",NULL);
243: if (nep->ops->setfromoptions) {
244: (*nep->ops->setfromoptions)(PetscOptionsObject,nep);
245: }
246: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)nep);
247: PetscOptionsEnd();
249: if (!nep->V) { NEPGetBV(nep,&nep->V); }
250: BVSetFromOptions(nep->V);
251: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
252: RGSetFromOptions(nep->rg);
253: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
254: DSSetFromOptions(nep->ds);
255: if (!nep->refineksp) { NEPRefineGetKSP(nep,&nep->refineksp); }
256: KSPSetFromOptions(nep->refineksp);
257: if (nep->fui==NEP_USER_INTERFACE_SPLIT) for (i=0;i<nep->nt;i++) {FNSetFromOptions(nep->f[i]);}
258: return(0);
259: }
261: /*@C
262: NEPGetTolerances - Gets the tolerance and maximum iteration count used
263: by the NEP convergence tests.
265: Not Collective
267: Input Parameter:
268: . nep - the nonlinear eigensolver context
270: Output Parameters:
271: + tol - the convergence tolerance
272: - maxits - maximum number of iterations
274: Notes:
275: The user can specify NULL for any parameter that is not needed.
277: Level: intermediate
279: .seealso: NEPSetTolerances()
280: @*/
281: PetscErrorCode NEPGetTolerances(NEP nep,PetscReal *tol,PetscInt *maxits)282: {
285: if (tol) *tol = nep->tol;
286: if (maxits) *maxits = nep->max_it;
287: return(0);
288: }
290: /*@
291: NEPSetTolerances - Sets the tolerance and maximum iteration count used
292: by the NEP convergence tests.
294: Logically Collective on NEP296: Input Parameters:
297: + nep - the nonlinear eigensolver context
298: . tol - the convergence tolerance
299: - maxits - maximum number of iterations to use
301: Options Database Keys:
302: + -nep_tol <tol> - Sets the convergence tolerance
303: - -nep_max_it <maxits> - Sets the maximum number of iterations allowed
305: Notes:
306: Use PETSC_DEFAULT for either argument to assign a reasonably good value.
308: Level: intermediate
310: .seealso: NEPGetTolerances()
311: @*/
312: PetscErrorCode NEPSetTolerances(NEP nep,PetscReal tol,PetscInt maxits)313: {
318: if (tol == PETSC_DEFAULT) {
319: nep->tol = PETSC_DEFAULT;
320: nep->state = NEP_STATE_INITIAL;
321: } else {
322: if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
323: nep->tol = tol;
324: }
325: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
326: nep->max_it = 0;
327: nep->state = NEP_STATE_INITIAL;
328: } else {
329: if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
330: nep->max_it = maxits;
331: }
332: return(0);
333: }
335: /*@C
336: NEPGetDimensions - Gets the number of eigenvalues to compute
337: and the dimension of the subspace.
339: Not Collective
341: Input Parameter:
342: . nep - the nonlinear eigensolver context
344: Output Parameters:
345: + nev - number of eigenvalues to compute
346: . ncv - the maximum dimension of the subspace to be used by the solver
347: - mpd - the maximum dimension allowed for the projected problem
349: Notes:
350: The user can specify NULL for any parameter that is not needed.
352: Level: intermediate
354: .seealso: NEPSetDimensions()
355: @*/
356: PetscErrorCode NEPGetDimensions(NEP nep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)357: {
360: if (nev) *nev = nep->nev;
361: if (ncv) *ncv = nep->ncv;
362: if (mpd) *mpd = nep->mpd;
363: return(0);
364: }
366: /*@
367: NEPSetDimensions - Sets the number of eigenvalues to compute
368: and the dimension of the subspace.
370: Logically Collective on NEP372: Input Parameters:
373: + nep - the nonlinear eigensolver context
374: . nev - number of eigenvalues to compute
375: . ncv - the maximum dimension of the subspace to be used by the solver
376: - mpd - the maximum dimension allowed for the projected problem
378: Options Database Keys:
379: + -nep_nev <nev> - Sets the number of eigenvalues
380: . -nep_ncv <ncv> - Sets the dimension of the subspace
381: - -nep_mpd <mpd> - Sets the maximum projected dimension
383: Notes:
384: Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
385: dependent on the solution method.
387: The parameters ncv and mpd are intimately related, so that the user is advised
388: to set one of them at most. Normal usage is that
389: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
390: (b) in cases where nev is large, the user sets mpd.
392: The value of ncv should always be between nev and (nev+mpd), typically
393: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
394: a smaller value should be used.
396: Level: intermediate
398: .seealso: NEPGetDimensions()
399: @*/
400: PetscErrorCode NEPSetDimensions(NEP nep,PetscInt nev,PetscInt ncv,PetscInt mpd)401: {
407: if (nev<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
408: nep->nev = nev;
409: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
410: nep->ncv = 0;
411: } else {
412: if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
413: nep->ncv = ncv;
414: }
415: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
416: nep->mpd = 0;
417: } else {
418: if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
419: nep->mpd = mpd;
420: }
421: nep->state = NEP_STATE_INITIAL;
422: return(0);
423: }
425: /*@
426: NEPSetWhichEigenpairs - Specifies which portion of the spectrum is
427: to be sought.
429: Logically Collective on NEP431: Input Parameters:
432: + nep - eigensolver context obtained from NEPCreate()
433: - which - the portion of the spectrum to be sought
435: Possible values:
436: The parameter 'which' can have one of these values
438: + NEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
439: . NEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
440: . NEP_LARGEST_REAL - largest real parts
441: . NEP_SMALLEST_REAL - smallest real parts
442: . NEP_LARGEST_IMAGINARY - largest imaginary parts
443: . NEP_SMALLEST_IMAGINARY - smallest imaginary parts
444: . NEP_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
445: . NEP_TARGET_REAL - eigenvalues with real part closest to target
446: . NEP_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
447: . NEP_ALL - all eigenvalues contained in a given region
448: - NEP_WHICH_USER - user defined ordering set with NEPSetEigenvalueComparison()
450: Options Database Keys:
451: + -nep_largest_magnitude - Sets largest eigenvalues in magnitude
452: . -nep_smallest_magnitude - Sets smallest eigenvalues in magnitude
453: . -nep_largest_real - Sets largest real parts
454: . -nep_smallest_real - Sets smallest real parts
455: . -nep_largest_imaginary - Sets largest imaginary parts
456: . -nep_smallest_imaginary - Sets smallest imaginary parts
457: . -nep_target_magnitude - Sets eigenvalues closest to target
458: . -nep_target_real - Sets real parts closest to target
459: . -nep_target_imaginary - Sets imaginary parts closest to target
460: - -nep_all - Sets all eigenvalues in a region
462: Notes:
463: Not all eigensolvers implemented in NEP account for all the possible values
464: stated above. If SLEPc is compiled for real numbers NEP_LARGEST_IMAGINARY465: and NEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
466: for eigenvalue selection.
468: The target is a scalar value provided with NEPSetTarget().
470: NEP_ALL is intended for use in the context of the CISS solver for
471: computing all eigenvalues in a region.
473: Level: intermediate
475: .seealso: NEPGetWhichEigenpairs(), NEPSetTarget(), NEPSetEigenvalueComparison(), NEPWhich476: @*/
477: PetscErrorCode NEPSetWhichEigenpairs(NEP nep,NEPWhich which)478: {
482: switch (which) {
483: case NEP_LARGEST_MAGNITUDE:
484: case NEP_SMALLEST_MAGNITUDE:
485: case NEP_LARGEST_REAL:
486: case NEP_SMALLEST_REAL:
487: case NEP_LARGEST_IMAGINARY:
488: case NEP_SMALLEST_IMAGINARY:
489: case NEP_TARGET_MAGNITUDE:
490: case NEP_TARGET_REAL:
491: #if defined(PETSC_USE_COMPLEX)
492: case NEP_TARGET_IMAGINARY:
493: #endif
494: case EPS_ALL:
495: case NEP_WHICH_USER:
496: if (nep->which != which) {
497: nep->state = NEP_STATE_INITIAL;
498: nep->which = which;
499: }
500: break;
501: default:502: SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
503: }
504: return(0);
505: }
507: /*@
508: NEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
509: sought.
511: Not Collective
513: Input Parameter:
514: . nep - eigensolver context obtained from NEPCreate()
516: Output Parameter:
517: . which - the portion of the spectrum to be sought
519: Notes:
520: See NEPSetWhichEigenpairs() for possible values of 'which'.
522: Level: intermediate
524: .seealso: NEPSetWhichEigenpairs(), NEPWhich525: @*/
526: PetscErrorCode NEPGetWhichEigenpairs(NEP nep,NEPWhich *which)527: {
531: *which = nep->which;
532: return(0);
533: }
535: /*@C
536: NEPSetEigenvalueComparison - Specifies the eigenvalue comparison function
537: when NEPSetWhichEigenpairs() is set to NEP_WHICH_USER.
539: Logically Collective on NEP541: Input Parameters:
542: + nep - eigensolver context obtained from NEPCreate()
543: . func - a pointer to the comparison function
544: - ctx - a context pointer (the last parameter to the comparison function)
546: Calling Sequence of func:
547: $ func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
549: + ar - real part of the 1st eigenvalue
550: . ai - imaginary part of the 1st eigenvalue
551: . br - real part of the 2nd eigenvalue
552: . bi - imaginary part of the 2nd eigenvalue
553: . res - result of comparison
554: - ctx - optional context, as set by NEPSetEigenvalueComparison()
556: Note:
557: The returning parameter 'res' can be
558: + negative - if the 1st eigenvalue is preferred to the 2st one
559: . zero - if both eigenvalues are equally preferred
560: - positive - if the 2st eigenvalue is preferred to the 1st one
562: Level: advanced
564: .seealso: NEPSetWhichEigenpairs(), NEPWhich565: @*/
566: PetscErrorCode NEPSetEigenvalueComparison(NEP nep,PetscErrorCode (*func)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)567: {
570: nep->sc->comparison = func;
571: nep->sc->comparisonctx = ctx;
572: nep->which = NEP_WHICH_USER;
573: return(0);
574: }
576: /*@
577: NEPSetProblemType - Specifies the type of the nonlinear eigenvalue problem.
579: Logically Collective on NEP581: Input Parameters:
582: + nep - the nonlinear eigensolver context
583: - type - a known type of nonlinear eigenvalue problem
585: Options Database Keys:
586: + -nep_general - general problem with no particular structure
587: - -nep_rational - a rational eigenvalue problem defined in split form with all f_i rational
589: Notes:
590: Allowed values for the problem type are: general (NEP_GENERAL), and rational
591: (NEP_RATIONAL).
593: This function is used to provide a hint to the NEP solver to exploit certain
594: properties of the nonlinear eigenproblem. This hint may be used or not,
595: depending on the solver. By default, no particular structure is assumed.
597: Level: intermediate
599: .seealso: NEPSetType(), NEPGetProblemType(), NEPProblemType600: @*/
601: PetscErrorCode NEPSetProblemType(NEP nep,NEPProblemType type)602: {
606: if (type!=NEP_GENERAL && type!=NEP_RATIONAL) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
607: if (type != nep->problem_type) {
608: nep->problem_type = type;
609: nep->state = NEP_STATE_INITIAL;
610: }
611: return(0);
612: }
614: /*@
615: NEPGetProblemType - Gets the problem type from the NEP object.
617: Not Collective
619: Input Parameter:
620: . nep - the nonlinear eigensolver context
622: Output Parameter:
623: . type - the problem type
625: Level: intermediate
627: .seealso: NEPSetProblemType(), NEPProblemType628: @*/
629: PetscErrorCode NEPGetProblemType(NEP nep,NEPProblemType *type)630: {
634: *type = nep->problem_type;
635: return(0);
636: }
638: /*@C
639: NEPSetConvergenceTestFunction - Sets a function to compute the error estimate
640: used in the convergence test.
642: Logically Collective on NEP644: Input Parameters:
645: + nep - nonlinear eigensolver context obtained from NEPCreate()
646: . func - a pointer to the convergence test function
647: . ctx - context for private data for the convergence routine (may be null)
648: - destroy - a routine for destroying the context (may be null)
650: Calling Sequence of func:
651: $ func(NEP nep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
653: + nep - nonlinear eigensolver context obtained from NEPCreate()
654: . eigr - real part of the eigenvalue
655: . eigi - imaginary part of the eigenvalue
656: . res - residual norm associated to the eigenpair
657: . errest - (output) computed error estimate
658: - ctx - optional context, as set by NEPSetConvergenceTestFunction()
660: Note:
661: If the error estimate returned by the convergence test function is less than
662: the tolerance, then the eigenvalue is accepted as converged.
664: Level: advanced
666: .seealso: NEPSetConvergenceTest(), NEPSetTolerances()
667: @*/
668: PetscErrorCode NEPSetConvergenceTestFunction(NEP nep,PetscErrorCode (*func)(NEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))669: {
674: if (nep->convergeddestroy) {
675: (*nep->convergeddestroy)(nep->convergedctx);
676: }
677: nep->convergeduser = func;
678: nep->convergeddestroy = destroy;
679: nep->convergedctx = ctx;
680: if (func == NEPConvergedRelative) nep->conv = NEP_CONV_REL;
681: else if (func == NEPConvergedNorm) nep->conv = NEP_CONV_NORM;
682: else if (func == NEPConvergedAbsolute) nep->conv = NEP_CONV_ABS;
683: else {
684: nep->conv = NEP_CONV_USER;
685: nep->converged = nep->convergeduser;
686: }
687: return(0);
688: }
690: /*@
691: NEPSetConvergenceTest - Specifies how to compute the error estimate
692: used in the convergence test.
694: Logically Collective on NEP696: Input Parameters:
697: + nep - nonlinear eigensolver context obtained from NEPCreate()
698: - conv - the type of convergence test
700: Options Database Keys:
701: + -nep_conv_abs - Sets the absolute convergence test
702: . -nep_conv_rel - Sets the convergence test relative to the eigenvalue
703: - -nep_conv_user - Selects the user-defined convergence test
705: Note:
706: The parameter 'conv' can have one of these values
707: + NEP_CONV_ABS - absolute error ||r||
708: . NEP_CONV_REL - error relative to the eigenvalue l, ||r||/|l|
709: . NEP_CONV_NORM - error relative matrix norms, ||r||/sum_i(|f_i(l)|*||A_i||)
710: - NEP_CONV_USER - function set by NEPSetConvergenceTestFunction()
712: Level: intermediate
714: .seealso: NEPGetConvergenceTest(), NEPSetConvergenceTestFunction(), NEPSetStoppingTest(), NEPConv715: @*/
716: PetscErrorCode NEPSetConvergenceTest(NEP nep,NEPConv conv)717: {
721: switch (conv) {
722: case NEP_CONV_ABS: nep->converged = NEPConvergedAbsolute; break;
723: case NEP_CONV_REL: nep->converged = NEPConvergedRelative; break;
724: case NEP_CONV_NORM: nep->converged = NEPConvergedNorm; break;
725: case NEP_CONV_USER:
726: if (!nep->convergeduser) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ORDER,"Must call NEPSetConvergenceTestFunction() first");
727: nep->converged = nep->convergeduser;
728: break;
729: default:730: SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
731: }
732: nep->conv = conv;
733: return(0);
734: }
736: /*@
737: NEPGetConvergenceTest - Gets the method used to compute the error estimate
738: used in the convergence test.
740: Not Collective
742: Input Parameters:
743: . nep - nonlinear eigensolver context obtained from NEPCreate()
745: Output Parameters:
746: . conv - the type of convergence test
748: Level: intermediate
750: .seealso: NEPSetConvergenceTest(), NEPConv751: @*/
752: PetscErrorCode NEPGetConvergenceTest(NEP nep,NEPConv *conv)753: {
757: *conv = nep->conv;
758: return(0);
759: }
761: /*@C
762: NEPSetStoppingTestFunction - Sets a function to decide when to stop the outer
763: iteration of the eigensolver.
765: Logically Collective on NEP767: Input Parameters:
768: + nep - nonlinear eigensolver context obtained from NEPCreate()
769: . func - pointer to the stopping test function
770: . ctx - context for private data for the stopping routine (may be null)
771: - destroy - a routine for destroying the context (may be null)
773: Calling Sequence of func:
774: $ func(NEP nep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,NEPConvergedReason *reason,void *ctx)
776: + nep - nonlinear eigensolver context obtained from NEPCreate()
777: . its - current number of iterations
778: . max_it - maximum number of iterations
779: . nconv - number of currently converged eigenpairs
780: . nev - number of requested eigenpairs
781: . reason - (output) result of the stopping test
782: - ctx - optional context, as set by NEPSetStoppingTestFunction()
784: Note:
785: Normal usage is to first call the default routine NEPStoppingBasic() and then
786: set reason to NEP_CONVERGED_USER if some user-defined conditions have been
787: met. To let the eigensolver continue iterating, the result must be left as
788: NEP_CONVERGED_ITERATING.
790: Level: advanced
792: .seealso: NEPSetStoppingTest(), NEPStoppingBasic()
793: @*/
794: PetscErrorCode NEPSetStoppingTestFunction(NEP nep,PetscErrorCode (*func)(NEP,PetscInt,PetscInt,PetscInt,PetscInt,NEPConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))795: {
800: if (nep->stoppingdestroy) {
801: (*nep->stoppingdestroy)(nep->stoppingctx);
802: }
803: nep->stoppinguser = func;
804: nep->stoppingdestroy = destroy;
805: nep->stoppingctx = ctx;
806: if (func == NEPStoppingBasic) nep->stop = NEP_STOP_BASIC;
807: else {
808: nep->stop = NEP_STOP_USER;
809: nep->stopping = nep->stoppinguser;
810: }
811: return(0);
812: }
814: /*@
815: NEPSetStoppingTest - Specifies how to decide the termination of the outer
816: loop of the eigensolver.
818: Logically Collective on NEP820: Input Parameters:
821: + nep - nonlinear eigensolver context obtained from NEPCreate()
822: - stop - the type of stopping test
824: Options Database Keys:
825: + -nep_stop_basic - Sets the default stopping test
826: - -nep_stop_user - Selects the user-defined stopping test
828: Note:
829: The parameter 'stop' can have one of these values
830: + NEP_STOP_BASIC - default stopping test
831: - NEP_STOP_USER - function set by NEPSetStoppingTestFunction()
833: Level: advanced
835: .seealso: NEPGetStoppingTest(), NEPSetStoppingTestFunction(), NEPSetConvergenceTest(), NEPStop836: @*/
837: PetscErrorCode NEPSetStoppingTest(NEP nep,NEPStop stop)838: {
842: switch (stop) {
843: case NEP_STOP_BASIC: nep->stopping = NEPStoppingBasic; break;
844: case NEP_STOP_USER:
845: if (!nep->stoppinguser) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ORDER,"Must call NEPSetStoppingTestFunction() first");
846: nep->stopping = nep->stoppinguser;
847: break;
848: default:849: SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
850: }
851: nep->stop = stop;
852: return(0);
853: }
855: /*@
856: NEPGetStoppingTest - Gets the method used to decide the termination of the outer
857: loop of the eigensolver.
859: Not Collective
861: Input Parameters:
862: . nep - nonlinear eigensolver context obtained from NEPCreate()
864: Output Parameters:
865: . stop - the type of stopping test
867: Level: advanced
869: .seealso: NEPSetStoppingTest(), NEPStop870: @*/
871: PetscErrorCode NEPGetStoppingTest(NEP nep,NEPStop *stop)872: {
876: *stop = nep->stop;
877: return(0);
878: }
880: /*@
881: NEPSetTrackAll - Specifies if the solver must compute the residual of all
882: approximate eigenpairs or not.
884: Logically Collective on NEP886: Input Parameters:
887: + nep - the eigensolver context
888: - trackall - whether compute all residuals or not
890: Notes:
891: If the user sets trackall=PETSC_TRUE then the solver explicitly computes
892: the residual for each eigenpair approximation. Computing the residual is
893: usually an expensive operation and solvers commonly compute the associated
894: residual to the first unconverged eigenpair.
895: The options '-nep_monitor_all' and '-nep_monitor_lg_all' automatically
896: activate this option.
898: Level: developer
900: .seealso: NEPGetTrackAll()
901: @*/
902: PetscErrorCode NEPSetTrackAll(NEP nep,PetscBool trackall)903: {
907: nep->trackall = trackall;
908: return(0);
909: }
911: /*@
912: NEPGetTrackAll - Returns the flag indicating whether all residual norms must
913: be computed or not.
915: Not Collective
917: Input Parameter:
918: . nep - the eigensolver context
920: Output Parameter:
921: . trackall - the returned flag
923: Level: developer
925: .seealso: NEPSetTrackAll()
926: @*/
927: PetscErrorCode NEPGetTrackAll(NEP nep,PetscBool *trackall)928: {
932: *trackall = nep->trackall;
933: return(0);
934: }
936: /*@
937: NEPSetRefine - Specifies the refinement type (and options) to be used
938: after the solve.
940: Logically Collective on NEP942: Input Parameters:
943: + nep - the nonlinear eigensolver context
944: . refine - refinement type
945: . npart - number of partitions of the communicator
946: . tol - the convergence tolerance
947: . its - maximum number of refinement iterations
948: - scheme - which scheme to be used for solving the involved linear systems
950: Options Database Keys:
951: + -nep_refine <type> - refinement type, one of <none,simple,multiple>
952: . -nep_refine_partitions <n> - the number of partitions
953: . -nep_refine_tol <tol> - the tolerance
954: . -nep_refine_its <its> - number of iterations
955: - -nep_refine_scheme - to set the scheme for the linear solves
957: Notes:
958: By default, iterative refinement is disabled, since it may be very
959: costly. There are two possible refinement strategies: simple and multiple.
960: The simple approach performs iterative refinement on each of the
961: converged eigenpairs individually, whereas the multiple strategy works
962: with the invariant pair as a whole, refining all eigenpairs simultaneously.
963: The latter may be required for the case of multiple eigenvalues.
965: In some cases, especially when using direct solvers within the
966: iterative refinement method, it may be helpful for improved scalability
967: to split the communicator in several partitions. The npart parameter
968: indicates how many partitions to use (defaults to 1).
970: The tol and its parameters specify the stopping criterion. In the simple
971: method, refinement continues until the residual of each eigenpair is
972: below the tolerance (tol defaults to the NEP tol, but may be set to a
973: different value). In contrast, the multiple method simply performs its
974: refinement iterations (just one by default).
976: The scheme argument is used to change the way in which linear systems are
977: solved. Possible choices are: explicit, mixed block elimination (MBE),
978: and Schur complement.
980: Level: intermediate
982: .seealso: NEPGetRefine()
983: @*/
984: PetscErrorCode NEPSetRefine(NEP nep,NEPRefine refine,PetscInt npart,PetscReal tol,PetscInt its,NEPRefineScheme scheme)985: {
987: PetscMPIInt size;
996: nep->refine = refine;
997: if (refine) { /* process parameters only if not REFINE_NONE */
998: if (npart!=nep->npart) {
999: PetscSubcommDestroy(&nep->refinesubc);
1000: KSPDestroy(&nep->refineksp);
1001: }
1002: if (npart == PETSC_DEFAULT || npart == PETSC_DECIDE) {
1003: nep->npart = 1;
1004: } else {
1005: MPI_Comm_size(PetscObjectComm((PetscObject)nep),&size);
1006: if (npart<1 || npart>size) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of npart");
1007: nep->npart = npart;
1008: }
1009: if (tol == PETSC_DEFAULT || tol == PETSC_DECIDE) {
1010: nep->rtol = PETSC_DEFAULT;
1011: } else {
1012: if (tol<=0.0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
1013: nep->rtol = tol;
1014: }
1015: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) {
1016: nep->rits = PETSC_DEFAULT;
1017: } else {
1018: if (its<0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be >= 0");
1019: nep->rits = its;
1020: }
1021: nep->scheme = scheme;
1022: }
1023: nep->state = NEP_STATE_INITIAL;
1024: return(0);
1025: }
1027: /*@C
1028: NEPGetRefine - Gets the refinement strategy used by the NEP object, and the
1029: associated parameters.
1031: Not Collective
1033: Input Parameter:
1034: . nep - the nonlinear eigensolver context
1036: Output Parameters:
1037: + refine - refinement type
1038: . npart - number of partitions of the communicator
1039: . tol - the convergence tolerance
1040: - its - maximum number of refinement iterations
1041: - scheme - the scheme used for solving linear systems
1043: Level: intermediate
1045: Note:
1046: The user can specify NULL for any parameter that is not needed.
1048: .seealso: NEPSetRefine()
1049: @*/
1050: PetscErrorCode NEPGetRefine(NEP nep,NEPRefine *refine,PetscInt *npart,PetscReal *tol,PetscInt *its,NEPRefineScheme *scheme)1051: {
1054: if (refine) *refine = nep->refine;
1055: if (npart) *npart = nep->npart;
1056: if (tol) *tol = nep->rtol;
1057: if (its) *its = nep->rits;
1058: if (scheme) *scheme = nep->scheme;
1059: return(0);
1060: }
1062: /*@C
1063: NEPSetOptionsPrefix - Sets the prefix used for searching for all
1064: NEP options in the database.
1066: Logically Collective on NEP1068: Input Parameters:
1069: + nep - the nonlinear eigensolver context
1070: - prefix - the prefix string to prepend to all NEP option requests
1072: Notes:
1073: A hyphen (-) must NOT be given at the beginning of the prefix name.
1074: The first character of all runtime options is AUTOMATICALLY the
1075: hyphen.
1077: For example, to distinguish between the runtime options for two
1078: different NEP contexts, one could call
1079: .vb
1080: NEPSetOptionsPrefix(nep1,"neig1_")
1081: NEPSetOptionsPrefix(nep2,"neig2_")
1082: .ve
1084: Level: advanced
1086: .seealso: NEPAppendOptionsPrefix(), NEPGetOptionsPrefix()
1087: @*/
1088: PetscErrorCode NEPSetOptionsPrefix(NEP nep,const char *prefix)1089: {
1094: if (!nep->V) { NEPGetBV(nep,&nep->V); }
1095: BVSetOptionsPrefix(nep->V,prefix);
1096: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
1097: DSSetOptionsPrefix(nep->ds,prefix);
1098: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
1099: RGSetOptionsPrefix(nep->rg,prefix);
1100: PetscObjectSetOptionsPrefix((PetscObject)nep,prefix);
1101: return(0);
1102: }
1104: /*@C
1105: NEPAppendOptionsPrefix - Appends to the prefix used for searching for all
1106: NEP options in the database.
1108: Logically Collective on NEP1110: Input Parameters:
1111: + nep - the nonlinear eigensolver context
1112: - prefix - the prefix string to prepend to all NEP option requests
1114: Notes:
1115: A hyphen (-) must NOT be given at the beginning of the prefix name.
1116: The first character of all runtime options is AUTOMATICALLY the hyphen.
1118: Level: advanced
1120: .seealso: NEPSetOptionsPrefix(), NEPGetOptionsPrefix()
1121: @*/
1122: PetscErrorCode NEPAppendOptionsPrefix(NEP nep,const char *prefix)1123: {
1128: if (!nep->V) { NEPGetBV(nep,&nep->V); }
1129: BVAppendOptionsPrefix(nep->V,prefix);
1130: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
1131: DSAppendOptionsPrefix(nep->ds,prefix);
1132: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
1133: RGAppendOptionsPrefix(nep->rg,prefix);
1134: PetscObjectAppendOptionsPrefix((PetscObject)nep,prefix);
1135: return(0);
1136: }
1138: /*@C
1139: NEPGetOptionsPrefix - Gets the prefix used for searching for all
1140: NEP options in the database.
1142: Not Collective
1144: Input Parameters:
1145: . nep - the nonlinear eigensolver context
1147: Output Parameters:
1148: . prefix - pointer to the prefix string used is returned
1150: Note:
1151: On the Fortran side, the user should pass in a string 'prefix' of
1152: sufficient length to hold the prefix.
1154: Level: advanced
1156: .seealso: NEPSetOptionsPrefix(), NEPAppendOptionsPrefix()
1157: @*/
1158: PetscErrorCode NEPGetOptionsPrefix(NEP nep,const char *prefix[])1159: {
1165: PetscObjectGetOptionsPrefix((PetscObject)nep,prefix);
1166: return(0);
1167: }