Actual source code: nepopts.c
1: /*
2: NEP routines related to options that can be set via the command-line
3: or procedurally.
5: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
6: SLEPc - Scalable Library for Eigenvalue Problem Computations
7: Copyright (c) 2002-2013, Universitat Politecnica de Valencia, Spain
9: This file is part of SLEPc.
11: SLEPc is free software: you can redistribute it and/or modify it under the
12: terms of version 3 of the GNU Lesser General Public License as published by
13: the Free Software Foundation.
15: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
16: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
18: more details.
20: You should have received a copy of the GNU Lesser General Public License
21: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
22: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
23: */
25: #include <slepc-private/nepimpl.h> /*I "slepcnep.h" I*/
29: /*@
30: NEPSetFromOptions - Sets NEP options from the options database.
31: This routine must be called before NEPSetUp() if the user is to be
32: allowed to set the solver type.
34: Collective on NEP
36: Input Parameters:
37: . nep - the nonlinear eigensolver context
39: Notes:
40: To see all options, run your program with the -help option.
42: Level: beginner
43: @*/
44: PetscErrorCode NEPSetFromOptions(NEP nep)
45: {
46: PetscErrorCode ierr;
47: char type[256],monfilename[PETSC_MAX_PATH_LEN];
48: PetscBool flg;
49: PetscReal r1,r2,r3;
50: PetscScalar s;
51: PetscInt i,j,k;
52: PetscViewer monviewer;
53: SlepcConvMonitor ctx;
54: MatStructure matflag;
58: if (!NEPRegisterAllCalled) { NEPRegisterAll(); }
59: PetscObjectOptionsBegin((PetscObject)nep);
60: PetscOptionsList("-nep_type","Nonlinear Eigenvalue Problem method","NEPSetType",NEPList,(char*)(((PetscObject)nep)->type_name?((PetscObject)nep)->type_name:NEPRII),type,256,&flg);
61: if (flg) {
62: NEPSetType(nep,type);
63: } else if (!((PetscObject)nep)->type_name) {
64: NEPSetType(nep,NEPRII);
65: }
67: r1 = r2 = r3 = i = j = 0;
68: PetscOptionsInt("-nep_max_it","Maximum number of iterations","NEPSetTolerances",nep->max_it,&i,NULL);
69: PetscOptionsInt("-nep_max_funcs","Maximum number of function evaluations","NEPSetTolerances",nep->max_funcs,&j,NULL);
70: PetscOptionsReal("-nep_atol","Absolute tolerance for residual norm","NEPSetTolerances",nep->abstol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:nep->abstol,&r1,NULL);
71: PetscOptionsReal("-nep_rtol","Relative tolerance for residual norm","NEPSetTolerances",nep->rtol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:nep->rtol,&r2,NULL);
72: PetscOptionsReal("-nep_stol","Relative tolerance for step length","NEPSetTolerances",nep->stol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:nep->stol,&r3,NULL);
73: NEPSetTolerances(nep,r1,r2,r3,i,j);
74: flg = PETSC_FALSE;
75: PetscOptionsBool("-nep_convergence_default","Default (relative error) convergence test","NEPSetConvergenceTest",flg,&flg,NULL);
76: if (flg) { NEPSetConvergenceTest(nep,NEPConvergedDefault,NULL,NULL); }
78: i = j = k = 0;
79: PetscOptionsInt("-nep_nev","Number of eigenvalues to compute","NEPSetDimensions",nep->nev,&i,NULL);
80: PetscOptionsInt("-nep_ncv","Number of basis vectors","NEPSetDimensions",nep->ncv,&j,NULL);
81: PetscOptionsInt("-nep_mpd","Maximum dimension of projected problem","NEPSetDimensions",nep->mpd,&k,NULL);
82: NEPSetDimensions(nep,i,j,k);
84: i = 0;
85: PetscOptionsInt("-nep_lag_preconditioner","Interval to rebuild preconditioner","NEPSetLagPreconditioner",nep->lag,&i,&flg);
86: if (flg) { NEPSetLagPreconditioner(nep,i); }
88: PetscOptionsBool("-nep_const_correction_tol","Constant correction tolerance for the linear solver","NEPSetConstCorrectionTol",nep->cctol,&nep->cctol,NULL);
90: PetscOptionsScalar("-nep_target","Value of the target","NEPSetTarget",nep->target,&s,&flg);
91: if (flg) {
92: NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE);
93: NEPSetTarget(nep,s);
94: }
96: /* -----------------------------------------------------------------------*/
97: /*
98: Cancels all monitors hardwired into code before call to NEPSetFromOptions()
99: */
100: flg = PETSC_FALSE;
101: PetscOptionsBool("-nep_monitor_cancel","Remove any hardwired monitor routines","NEPMonitorCancel",flg,&flg,NULL);
102: if (flg) {
103: NEPMonitorCancel(nep);
104: }
105: /*
106: Prints approximate eigenvalues and error estimates at each iteration
107: */
108: PetscOptionsString("-nep_monitor","Monitor first unconverged approximate eigenvalue and error estimate","NEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
109: if (flg) {
110: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)nep),monfilename,&monviewer);
111: NEPMonitorSet(nep,NEPMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
112: }
113: PetscOptionsString("-nep_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","NEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
114: if (flg) {
115: PetscNew(struct _n_SlepcConvMonitor,&ctx);
116: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)nep),monfilename,&ctx->viewer);
117: NEPMonitorSet(nep,NEPMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
118: }
119: PetscOptionsString("-nep_monitor_all","Monitor approximate eigenvalues and error estimates","NEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
120: if (flg) {
121: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)nep),monfilename,&monviewer);
122: NEPMonitorSet(nep,NEPMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
123: NEPSetTrackAll(nep,PETSC_TRUE);
124: }
125: flg = PETSC_FALSE;
126: PetscOptionsBool("-nep_monitor_lg","Monitor first unconverged approximate error estimate graphically","NEPMonitorSet",flg,&flg,NULL);
127: if (flg) {
128: NEPMonitorSet(nep,NEPMonitorLG,NULL,NULL);
129: }
130: flg = PETSC_FALSE;
131: PetscOptionsBool("-nep_monitor_lg_all","Monitor error estimates graphically","NEPMonitorSet",flg,&flg,NULL);
132: if (flg) {
133: NEPMonitorSet(nep,NEPMonitorLGAll,NULL,NULL);
134: NEPSetTrackAll(nep,PETSC_TRUE);
135: }
136: /* -----------------------------------------------------------------------*/
138: PetscOptionsBoolGroupBegin("-nep_largest_magnitude","compute largest eigenvalues in magnitude","NEPSetWhichEigenpairs",&flg);
139: if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_MAGNITUDE); }
140: PetscOptionsBoolGroup("-nep_smallest_magnitude","compute smallest eigenvalues in magnitude","NEPSetWhichEigenpairs",&flg);
141: if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_MAGNITUDE); }
142: PetscOptionsBoolGroup("-nep_largest_real","compute largest real parts","NEPSetWhichEigenpairs",&flg);
143: if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_REAL); }
144: PetscOptionsBoolGroup("-nep_smallest_real","compute smallest real parts","NEPSetWhichEigenpairs",&flg);
145: if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_REAL); }
146: PetscOptionsBoolGroup("-nep_largest_imaginary","compute largest imaginary parts","NEPSetWhichEigenpairs",&flg);
147: if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_IMAGINARY); }
148: PetscOptionsBoolGroup("-nep_smallest_imaginary","compute smallest imaginary parts","NEPSetWhichEigenpairs",&flg);
149: if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_IMAGINARY); }
150: PetscOptionsBoolGroup("-nep_target_magnitude","compute nearest eigenvalues to target","NEPSetWhichEigenpairs",&flg);
151: if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE); }
152: PetscOptionsBoolGroup("-nep_target_real","compute eigenvalues with real parts close to target","NEPSetWhichEigenpairs",&flg);
153: if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_REAL); }
154: PetscOptionsBoolGroupEnd("-nep_target_imaginary","compute eigenvalues with imaginary parts close to target","NEPSetWhichEigenpairs",&flg);
155: if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_IMAGINARY); }
157: PetscOptionsName("-nep_view","Print detailed information on solver used","NEPView",0);
158: PetscOptionsName("-nep_plot_eigs","Make a plot of the computed eigenvalues","NEPSolve",0);
160: if (nep->ops->setfromoptions) {
161: (*nep->ops->setfromoptions)(nep);
162: }
163: PetscObjectProcessOptionsHandlers((PetscObject)nep);
164: PetscOptionsEnd();
166: if (!nep->ip) { NEPGetIP(nep,&nep->ip); }
167: IPSetFromOptions(nep->ip);
168: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
169: DSSetFromOptions(nep->ds);
170: if (!nep->ksp) { NEPGetKSP(nep,&nep->ksp); }
171: KSPGetOperators(nep->ksp,NULL,NULL,&matflag);
172: KSPSetOperators(nep->ksp,nep->function,nep->function_pre,matflag);
173: KSPSetFromOptions(nep->ksp);
174: PetscRandomSetFromOptions(nep->rand);
175: return(0);
176: }
180: /*@
181: NEPGetTolerances - Gets the tolerance and maximum iteration count used
182: by the NEP convergence tests.
184: Not Collective
186: Input Parameter:
187: . nep - the nonlinear eigensolver context
189: Output Parameters:
190: + abstol - absolute convergence tolerance
191: . rtol - relative convergence tolerance
192: . stol - convergence tolerance in terms of the norm of the change in the
193: solution between steps, || delta x || < stol*|| x ||
194: . maxit - maximum number of iterations
195: - maxf - maximum number of function evaluations
197: Notes:
198: The user can specify NULL for any parameter that is not needed.
200: Level: intermediate
202: .seealso: NEPSetTolerances()
203: @*/
204: PetscErrorCode NEPGetTolerances(NEP nep,PetscReal *abstol,PetscReal *rtol,PetscReal *stol,PetscInt *maxit,PetscInt *maxf)
205: {
208: if (abstol) *abstol = nep->abstol;
209: if (rtol) *rtol = nep->rtol;
210: if (stol) *stol = nep->stol;
211: if (maxit) *maxit = nep->max_it;
212: if (maxf) *maxf = nep->max_funcs;
213: return(0);
214: }
218: /*@
219: NEPSetTolerances - Sets various parameters used in convergence tests.
221: Logically Collective on NEP
223: Input Parameters:
224: + nep - the nonlinear eigensolver context
225: . abstol - absolute convergence tolerance
226: . rtol - relative convergence tolerance
227: . stol - convergence tolerance in terms of the norm of the change in the
228: solution between steps, || delta x || < stol*|| x ||
229: . maxit - maximum number of iterations
230: - maxf - maximum number of function evaluations
232: Options Database Keys:
233: + -nep_atol <abstol> - Sets abstol
234: . -nep_rtol <rtol> - Sets rtol
235: . -nep_stol <stol> - Sets stol
236: . -nep_max_it <maxit> - Sets maxit
237: - -nep_max_funcs <maxf> - Sets maxf
239: Notes:
240: Pass 0 for an argument that need not be changed.
242: Use PETSC_DECIDE for maxits to assign a reasonably good value, which is
243: dependent on the solution method.
245: Level: intermediate
247: .seealso: NEPGetTolerances()
248: @*/
249: PetscErrorCode NEPSetTolerances(NEP nep,PetscReal abstol,PetscReal rtol,PetscReal stol,PetscInt maxit,PetscInt maxf)
250: {
258: if (abstol) {
259: if (abstol == PETSC_DEFAULT) {
260: nep->abstol = PETSC_DEFAULT;
261: } else {
262: if (abstol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %G must be non-negative",abstol);
263: nep->abstol = abstol;
264: }
265: }
266: if (rtol) {
267: if (rtol == PETSC_DEFAULT) {
268: nep->rtol = PETSC_DEFAULT;
269: } else {
270: if (rtol < 0.0 || 1.0 <= rtol) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Relative tolerance %G must be non-negative and less than 1.0",rtol);
271: nep->rtol = rtol;
272: }
273: }
274: if (stol) {
275: if (stol == PETSC_DEFAULT) {
276: nep->stol = PETSC_DEFAULT;
277: } else {
278: if (stol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Step tolerance %G must be non-negative",stol);
279: nep->stol = stol;
280: }
281: }
282: if (maxit) {
283: if (maxit == PETSC_DEFAULT || maxit == PETSC_DECIDE) {
284: nep->max_it = 0;
285: nep->setupcalled = 0;
286: } else {
287: if (maxit < 0) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",maxit);
288: nep->max_it = maxit;
289: }
290: }
291: if (maxf) {
292: if (maxf == PETSC_DEFAULT || maxf == PETSC_DECIDE) {
293: nep->max_it = 0;
294: nep->setupcalled = 0;
295: } else {
296: if (maxf < 0) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of function evaluations %D must be non-negative",maxf);
297: nep->max_funcs = maxf;
298: }
299: }
300: return(0);
301: }
305: /*@
306: NEPGetDimensions - Gets the number of eigenvalues to compute
307: and the dimension of the subspace.
309: Not Collective
311: Input Parameter:
312: . nep - the nonlinear eigensolver context
314: Output Parameters:
315: + nev - number of eigenvalues to compute
316: . ncv - the maximum dimension of the subspace to be used by the solver
317: - mpd - the maximum dimension allowed for the projected problem
319: Notes:
320: The user can specify NULL for any parameter that is not needed.
322: Level: intermediate
324: .seealso: NEPSetDimensions()
325: @*/
326: PetscErrorCode NEPGetDimensions(NEP nep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
327: {
330: if (nev) *nev = nep->nev;
331: if (ncv) *ncv = nep->ncv;
332: if (mpd) *mpd = nep->mpd;
333: return(0);
334: }
338: /*@
339: NEPSetDimensions - Sets the number of eigenvalues to compute
340: and the dimension of the subspace.
342: Logically Collective on NEP
344: Input Parameters:
345: + nep - the nonlinear eigensolver context
346: . nev - number of eigenvalues to compute
347: . ncv - the maximum dimension of the subspace to be used by the solver
348: - mpd - the maximum dimension allowed for the projected problem
350: Options Database Keys:
351: + -nep_nev <nev> - Sets the number of eigenvalues
352: . -nep_ncv <ncv> - Sets the dimension of the subspace
353: - -nep_mpd <mpd> - Sets the maximum projected dimension
355: Notes:
356: Pass 0 to retain the previous value of any parameter.
358: Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
359: dependent on the solution method.
361: The parameters ncv and mpd are intimately related, so that the user is advised
362: to set one of them at most. Normal usage is that
363: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
364: (b) in cases where nev is large, the user sets mpd.
366: The value of ncv should always be between nev and (nev+mpd), typically
367: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
368: a smaller value should be used.
370: Level: intermediate
372: .seealso: NEPGetDimensions()
373: @*/
374: PetscErrorCode NEPSetDimensions(NEP nep,PetscInt nev,PetscInt ncv,PetscInt mpd)
375: {
381: if (nev) {
382: if (nev<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
383: nep->nev = nev;
384: nep->setupcalled = 0;
385: }
386: if (ncv) {
387: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
388: nep->ncv = 0;
389: } else {
390: if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
391: nep->ncv = ncv;
392: }
393: nep->setupcalled = 0;
394: }
395: if (mpd) {
396: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
397: nep->mpd = 0;
398: } else {
399: if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
400: nep->mpd = mpd;
401: }
402: }
403: return(0);
404: }
408: /*@
409: NEPSetWhichEigenpairs - Specifies which portion of the spectrum is
410: to be sought.
412: Logically Collective on NEP
414: Input Parameters:
415: + nep - eigensolver context obtained from NEPCreate()
416: - which - the portion of the spectrum to be sought
418: Possible values:
419: The parameter 'which' can have one of these values
421: + NEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
422: . NEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
423: . NEP_LARGEST_REAL - largest real parts
424: . NEP_SMALLEST_REAL - smallest real parts
425: . NEP_LARGEST_IMAGINARY - largest imaginary parts
426: . NEP_SMALLEST_IMAGINARY - smallest imaginary parts
427: . NEP_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
428: . NEP_TARGET_REAL - eigenvalues with real part closest to target
429: - NEP_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
431: Options Database Keys:
432: + -nep_largest_magnitude - Sets largest eigenvalues in magnitude
433: . -nep_smallest_magnitude - Sets smallest eigenvalues in magnitude
434: . -nep_largest_real - Sets largest real parts
435: . -nep_smallest_real - Sets smallest real parts
436: . -nep_largest_imaginary - Sets largest imaginary parts
437: . -nep_smallest_imaginary - Sets smallest imaginary parts
438: . -nep_target_magnitude - Sets eigenvalues closest to target
439: . -nep_target_real - Sets real parts closest to target
440: - -nep_target_imaginary - Sets imaginary parts closest to target
442: Notes:
443: Not all eigensolvers implemented in NEP account for all the possible values
444: stated above. If SLEPc is compiled for real numbers NEP_LARGEST_IMAGINARY
445: and NEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
446: for eigenvalue selection.
448: Level: intermediate
450: .seealso: NEPGetWhichEigenpairs(), NEPWhich
451: @*/
452: PetscErrorCode NEPSetWhichEigenpairs(NEP nep,NEPWhich which)
453: {
457: if (which) {
458: if (which==PETSC_DECIDE || which==PETSC_DEFAULT) nep->which = (NEPWhich)0;
459: else switch (which) {
460: case NEP_LARGEST_MAGNITUDE:
461: case NEP_SMALLEST_MAGNITUDE:
462: case NEP_LARGEST_REAL:
463: case NEP_SMALLEST_REAL:
464: case NEP_LARGEST_IMAGINARY:
465: case NEP_SMALLEST_IMAGINARY:
466: case NEP_TARGET_MAGNITUDE:
467: case NEP_TARGET_REAL:
468: #if defined(PETSC_USE_COMPLEX)
469: case NEP_TARGET_IMAGINARY:
470: #endif
471: if (nep->which != which) {
472: nep->setupcalled = 0;
473: nep->which = which;
474: }
475: break;
476: default:
477: SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
478: }
479: }
480: return(0);
481: }
485: /*@C
486: NEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
487: sought.
489: Not Collective
491: Input Parameter:
492: . nep - eigensolver context obtained from NEPCreate()
494: Output Parameter:
495: . which - the portion of the spectrum to be sought
497: Notes:
498: See NEPSetWhichEigenpairs() for possible values of 'which'.
500: Level: intermediate
502: .seealso: NEPSetWhichEigenpairs(), NEPWhich
503: @*/
504: PetscErrorCode NEPGetWhichEigenpairs(NEP nep,NEPWhich *which)
505: {
509: *which = nep->which;
510: return(0);
511: }
515: /*@
516: NEPSetLagPreconditioner - Determines when the preconditioner is rebuilt in the
517: nonlinear solve.
519: Logically Collective on NEP
521: Input Parameters:
522: + nep - the NEP context
523: - lag - 0 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is
524: computed within the nonlinear iteration, 2 means every second time
525: the Jacobian is built, etc.
527: Options Database Keys:
528: . -nep_lag_preconditioner <lag>
530: Notes:
531: The default is 1.
532: The preconditioner is ALWAYS built in the first iteration of a nonlinear solve.
534: Level: intermediate
536: .seealso: NEPGetLagPreconditioner()
537: @*/
538: PetscErrorCode NEPSetLagPreconditioner(NEP nep,PetscInt lag)
539: {
543: if (lag<0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag must be non-negative");
544: nep->lag = lag;
545: return(0);
546: }
550: /*@
551: NEPGetLagPreconditioner - Indicates how often the preconditioner is rebuilt.
553: Not Collective
555: Input Parameter:
556: . nep - the NEP context
558: Output Parameter:
559: . lag - the lag parameter
561: Level: intermediate
563: .seealso: NEPSetLagPreconditioner()
564: @*/
565: PetscErrorCode NEPGetLagPreconditioner(NEP nep,PetscInt *lag)
566: {
570: *lag = nep->lag;
571: return(0);
572: }
576: /*@
577: NEPSetConstCorrectionTol - Sets a flag to keep the tolerance used
578: in the linear solver constant.
580: Logically Collective on NEP
582: Input Parameters:
583: + nep - the NEP context
584: - cct - a boolean value
586: Options Database Keys:
587: . -nep_const_correction_tol <cct>
589: Notes:
590: By default, an exponentially decreasing tolerance is set in the KSP used
591: within the nonlinear iteration, so that each Newton iteration requests
592: better accuracy than the previous one. The constant correction tolerance
593: flag stops this behaviour.
595: Level: intermediate
597: .seealso: NEPGetConstCorrectionTol()
598: @*/
599: PetscErrorCode NEPSetConstCorrectionTol(NEP nep,PetscBool cct)
600: {
604: nep->cctol = cct;
605: return(0);
606: }
610: /*@
611: NEPGetConstCorrectionTol - Returns the constant tolerance flag.
613: Not Collective
615: Input Parameter:
616: . nep - the NEP context
618: Output Parameter:
619: . cct - the value of the constant tolerance flag
621: Level: intermediate
623: .seealso: NEPSetConstCorrectionTol()
624: @*/
625: PetscErrorCode NEPGetConstCorrectionTol(NEP nep,PetscBool *cct)
626: {
630: *cct = nep->cctol;
631: return(0);
632: }
636: /*@C
637: NEPSetConvergenceTest - Sets the function to be used to test convergence
638: of the nonlinear iterative solution.
640: Logically Collective on NEP
642: Input Parameters:
643: + nep - the NEP context
644: . func - a pointer to the convergence test function
645: . ctx - [optional] context for private data for the convergence routine
646: (may be NULL)
647: - destroy - [optional] destructor for the context (may be NULL;
648: PETSC_NULL_FUNCTION in Fortran)
650: Calling Sequence of func:
651: $ func(NEP nep,PetscInt it,PetscReal xnorm,PetscReal snorm,PetscReal fnorm,NEPConvergedReason reason*,void *fctx)
653: + nep - the NEP context
654: . it - iteration number
655: . xnorm - norm of the current solution
656: . snorm - norm of the step (difference between two consecutive solutions)
657: . fnorm - norm of the function (residual)
658: . reason - (output) result of the convergence test
659: - fctx - optional context, as set by NEPSetConvergenceTest()
661: Level: advanced
663: .seealso: NEPSetTolerances()
664: @*/
665: PetscErrorCode NEPSetConvergenceTest(NEP nep,PetscErrorCode (*func)(NEP,PetscInt,PetscReal,PetscReal,PetscReal,NEPConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
666: {
671: if (nep->convergeddestroy) {
672: (*nep->convergeddestroy)(nep->convergedctx);
673: }
674: nep->converged = func;
675: nep->convergeddestroy = destroy;
676: nep->convergedctx = ctx;
677: return(0);
678: }
682: /*@
683: NEPSetTrackAll - Specifies if the solver must compute the residual of all
684: approximate eigenpairs or not.
686: Logically Collective on NEP
688: Input Parameters:
689: + nep - the eigensolver context
690: - trackall - whether compute all residuals or not
692: Notes:
693: If the user sets trackall=PETSC_TRUE then the solver explicitly computes
694: the residual for each eigenpair approximation. Computing the residual is
695: usually an expensive operation and solvers commonly compute the associated
696: residual to the first unconverged eigenpair.
697: The options '-nep_monitor_all' and '-nep_monitor_lg_all' automatically
698: activate this option.
700: Level: intermediate
702: .seealso: NEPGetTrackAll()
703: @*/
704: PetscErrorCode NEPSetTrackAll(NEP nep,PetscBool trackall)
705: {
709: nep->trackall = trackall;
710: return(0);
711: }
715: /*@
716: NEPGetTrackAll - Returns the flag indicating whether all residual norms must
717: be computed or not.
719: Not Collective
721: Input Parameter:
722: . nep - the eigensolver context
724: Output Parameter:
725: . trackall - the returned flag
727: Level: intermediate
729: .seealso: NEPSetTrackAll()
730: @*/
731: PetscErrorCode NEPGetTrackAll(NEP nep,PetscBool *trackall)
732: {
736: *trackall = nep->trackall;
737: return(0);
738: }
742: /*@C
743: NEPSetOptionsPrefix - Sets the prefix used for searching for all
744: NEP options in the database.
746: Logically Collective on NEP
748: Input Parameters:
749: + nep - the nonlinear eigensolver context
750: - prefix - the prefix string to prepend to all NEP option requests
752: Notes:
753: A hyphen (-) must NOT be given at the beginning of the prefix name.
754: The first character of all runtime options is AUTOMATICALLY the
755: hyphen.
757: For example, to distinguish between the runtime options for two
758: different NEP contexts, one could call
759: .vb
760: NEPSetOptionsPrefix(nep1,"neig1_")
761: NEPSetOptionsPrefix(nep2,"neig2_")
762: .ve
764: Level: advanced
766: .seealso: NEPAppendOptionsPrefix(), NEPGetOptionsPrefix()
767: @*/
768: PetscErrorCode NEPSetOptionsPrefix(NEP nep,const char *prefix)
769: {
774: if (!nep->ip) { NEPGetIP(nep,&nep->ip); }
775: IPSetOptionsPrefix(nep->ip,prefix);
776: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
777: DSSetOptionsPrefix(nep->ds,prefix);
778: if (!nep->ksp) { NEPGetKSP(nep,&nep->ksp); }
779: KSPSetOptionsPrefix(nep->ksp,prefix);
780: KSPAppendOptionsPrefix(nep->ksp,"nep_");
781: PetscObjectSetOptionsPrefix((PetscObject)nep,prefix);
782: return(0);
783: }
787: /*@C
788: NEPAppendOptionsPrefix - Appends to the prefix used for searching for all
789: NEP options in the database.
791: Logically Collective on NEP
793: Input Parameters:
794: + nep - the nonlinear eigensolver context
795: - prefix - the prefix string to prepend to all NEP option requests
797: Notes:
798: A hyphen (-) must NOT be given at the beginning of the prefix name.
799: The first character of all runtime options is AUTOMATICALLY the hyphen.
801: Level: advanced
803: .seealso: NEPSetOptionsPrefix(), NEPGetOptionsPrefix()
804: @*/
805: PetscErrorCode NEPAppendOptionsPrefix(NEP nep,const char *prefix)
806: {
811: if (!nep->ip) { NEPGetIP(nep,&nep->ip); }
812: IPSetOptionsPrefix(nep->ip,prefix);
813: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
814: DSSetOptionsPrefix(nep->ds,prefix);
815: if (!nep->ksp) { NEPGetKSP(nep,&nep->ksp); }
816: KSPSetOptionsPrefix(nep->ksp,prefix);
817: KSPAppendOptionsPrefix(nep->ksp,"nep_");
818: PetscObjectAppendOptionsPrefix((PetscObject)nep,prefix);
819: return(0);
820: }
824: /*@C
825: NEPGetOptionsPrefix - Gets the prefix used for searching for all
826: NEP options in the database.
828: Not Collective
830: Input Parameters:
831: . nep - the nonlinear eigensolver context
833: Output Parameters:
834: . prefix - pointer to the prefix string used is returned
836: Notes: On the fortran side, the user should pass in a string 'prefix' of
837: sufficient length to hold the prefix.
839: Level: advanced
841: .seealso: NEPSetOptionsPrefix(), NEPAppendOptionsPrefix()
842: @*/
843: PetscErrorCode NEPGetOptionsPrefix(NEP nep,const char *prefix[])
844: {
850: PetscObjectGetOptionsPrefix((PetscObject)nep,prefix);
851: return(0);
852: }