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-2016, 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*/
26: #include <petscdraw.h>
30: /*@C
31: NEPMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
32: indicated by the user.
34: Collective on NEP 36: Input Parameters:
37: + nep - the nonlinear eigensolver context
38: . name - the monitor option name
39: . help - message indicating what monitoring is done
40: . manual - manual page for the monitor
41: . monitor - the monitor function, whose context is a PetscViewerAndFormat
42: - trackall - whether this monitor tracks all eigenvalues or not
44: Level: developer
46: .seealso: NEPMonitorSet(), NEPSetTrackAll(), NEPConvMonitorSetFromOptions()
47: @*/
48: 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) 49: {
50: PetscErrorCode ierr;
51: PetscBool flg;
52: PetscViewer viewer;
53: PetscViewerFormat format;
54: PetscViewerAndFormat *vf;
57: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,name,&viewer,&format,&flg);
58: if (flg) {
59: PetscViewerAndFormatCreate(viewer,format,&vf);
60: PetscObjectDereference((PetscObject)viewer);
61: NEPMonitorSet(nep,(PetscErrorCode (*)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))monitor,vf,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
62: if (trackall) {
63: NEPSetTrackAll(nep,PETSC_TRUE);
64: }
65: }
66: return(0);
67: }
71: /*@C
72: NEPConvMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
73: indicated by the user (for monitors that only show iteration numbers of convergence).
75: Collective on NEP 77: Input Parameters:
78: + nep - the nonlinear eigensolver context
79: . name - the monitor option name
80: . help - message indicating what monitoring is done
81: . manual - manual page for the monitor
82: - monitor - the monitor function, whose context is a SlepcConvMonitor
84: Level: developer
86: .seealso: NEPMonitorSet(), NEPMonitorSetFromOptions()
87: @*/
88: PetscErrorCode NEPConvMonitorSetFromOptions(NEP nep,const char name[],const char help[],const char manual[],PetscErrorCode (*monitor)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,SlepcConvMonitor)) 89: {
90: PetscErrorCode ierr;
91: PetscBool flg;
92: PetscViewer viewer;
93: PetscViewerFormat format;
94: SlepcConvMonitor ctx;
97: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,name,&viewer,&format,&flg);
98: if (flg) {
99: SlepcConvMonitorCreate(viewer,format,&ctx);
100: PetscObjectDereference((PetscObject)viewer);
101: NEPMonitorSet(nep,(PetscErrorCode (*)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))monitor,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
102: }
103: return(0);
104: }
108: /*@
109: NEPSetFromOptions - Sets NEP options from the options database.
110: This routine must be called before NEPSetUp() if the user is to be
111: allowed to set the solver type.
113: Collective on NEP115: Input Parameters:
116: . nep - the nonlinear eigensolver context
118: Notes:
119: To see all options, run your program with the -help option.
121: Level: beginner
122: @*/
123: PetscErrorCode NEPSetFromOptions(NEP nep)124: {
126: char type[256];
127: PetscBool set,flg,flg1,flg2,flg3;
128: PetscReal r;
129: PetscScalar s;
130: PetscInt i,j,k;
131: PetscDrawLG lg;
135: NEPRegisterAll();
136: PetscObjectOptionsBegin((PetscObject)nep);
137: PetscOptionsFList("-nep_type","Nonlinear Eigenvalue Problem method","NEPSetType",NEPList,(char*)(((PetscObject)nep)->type_name?((PetscObject)nep)->type_name:NEPRII),type,256,&flg);
138: if (flg) {
139: NEPSetType(nep,type);
140: } else if (!((PetscObject)nep)->type_name) {
141: NEPSetType(nep,NEPRII);
142: }
144: PetscOptionsEnum("-nep_refine","Iterative refinement method","NEPSetRefine",NEPRefineTypes,(PetscEnum)nep->refine,(PetscEnum*)&nep->refine,NULL);
146: i = nep->npart;
147: PetscOptionsInt("-nep_refine_partitions","Number of partitions of the communicator for iterative refinement","NEPSetRefine",nep->npart,&i,&flg1);
148: r = nep->rtol;
149: PetscOptionsReal("-nep_refine_tol","Tolerance for iterative refinement","NEPSetRefine",nep->rtol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL/1000:nep->rtol,&r,&flg2);
150: j = nep->rits;
151: PetscOptionsInt("-nep_refine_its","Maximum number of iterations for iterative refinement","NEPSetRefine",nep->rits,&j,&flg3);
152: if (flg1 || flg2 || flg3) {
153: NEPSetRefine(nep,nep->refine,i,r,j,nep->scheme);
154: }
156: PetscOptionsEnum("-nep_refine_scheme","Scheme used for linear systems within iterative refinement","NEPSetRefine",NEPRefineSchemes,(PetscEnum)nep->scheme,(PetscEnum*)&nep->scheme,NULL);
158: i = nep->max_it? nep->max_it: PETSC_DEFAULT;
159: PetscOptionsInt("-nep_max_it","Maximum number of iterations","NEPSetTolerances",nep->max_it,&i,&flg1);
160: r = nep->tol;
161: PetscOptionsReal("-nep_tol","Tolerance","NEPSetTolerances",nep->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:nep->tol,&r,&flg2);
162: if (flg1 || flg2) {
163: NEPSetTolerances(nep,r,i);
164: }
166: PetscOptionsBoolGroupBegin("-nep_conv_rel","Relative error convergence test","NEPSetConvergenceTest",&flg);
167: if (flg) { NEPSetConvergenceTest(nep,NEP_CONV_REL); }
168: PetscOptionsBoolGroup("-nep_conv_norm","Convergence test relative to the matrix norms","NEPSetConvergenceTest",&flg);
169: if (flg) { NEPSetConvergenceTest(nep,NEP_CONV_NORM); }
170: PetscOptionsBoolGroup("-nep_conv_abs","Absolute error convergence test","NEPSetConvergenceTest",&flg);
171: if (flg) { NEPSetConvergenceTest(nep,NEP_CONV_ABS); }
172: PetscOptionsBoolGroupEnd("-nep_conv_user","User-defined convergence test","NEPSetConvergenceTest",&flg);
173: if (flg) { NEPSetConvergenceTest(nep,NEP_CONV_USER); }
175: PetscOptionsBoolGroupBegin("-nep_stop_basic","Stop iteration if all eigenvalues converged or max_it reached","NEPSetStoppingTest",&flg);
176: if (flg) { NEPSetStoppingTest(nep,NEP_STOP_BASIC); }
177: PetscOptionsBoolGroupEnd("-nep_stop_user","User-defined stopping test","NEPSetStoppingTest",&flg);
178: if (flg) { NEPSetStoppingTest(nep,NEP_STOP_USER); }
180: i = nep->nev;
181: PetscOptionsInt("-nep_nev","Number of eigenvalues to compute","NEPSetDimensions",nep->nev,&i,&flg1);
182: j = nep->ncv? nep->ncv: PETSC_DEFAULT;
183: PetscOptionsInt("-nep_ncv","Number of basis vectors","NEPSetDimensions",nep->ncv,&j,&flg2);
184: k = nep->mpd? nep->mpd: PETSC_DEFAULT;
185: PetscOptionsInt("-nep_mpd","Maximum dimension of projected problem","NEPSetDimensions",nep->mpd,&k,&flg3);
186: if (flg1 || flg2 || flg3) {
187: NEPSetDimensions(nep,i,j,k);
188: }
190: PetscOptionsScalar("-nep_target","Value of the target","NEPSetTarget",nep->target,&s,&flg);
191: if (flg) {
192: NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE);
193: NEPSetTarget(nep,s);
194: }
196: /* -----------------------------------------------------------------------*/
197: /*
198: Cancels all monitors hardwired into code before call to NEPSetFromOptions()
199: */
200: PetscOptionsBool("-nep_monitor_cancel","Remove any hardwired monitor routines","NEPMonitorCancel",PETSC_FALSE,&flg,&set);
201: if (set && flg) {
202: NEPMonitorCancel(nep);
203: }
204: /*
205: Text monitors
206: */
207: NEPMonitorSetFromOptions(nep,"-nep_monitor","Monitor first unconverged approximate eigenvalue and error estimate","NEPMonitorFirst",NEPMonitorFirst,PETSC_FALSE);
208: NEPConvMonitorSetFromOptions(nep,"-nep_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","NEPMonitorConverged",NEPMonitorConverged);
209: NEPMonitorSetFromOptions(nep,"-nep_monitor_all","Monitor approximate eigenvalues and error estimates","NEPMonitorAll",NEPMonitorAll,PETSC_TRUE);
210: /*
211: Line graph monitors
212: */
213: PetscOptionsBool("-nep_monitor_lg","Monitor first unconverged approximate error estimate graphically","NEPMonitorSet",PETSC_FALSE,&flg,&set);
214: if (set && flg) {
215: NEPMonitorLGCreate(PetscObjectComm((PetscObject)nep),NULL,"Error estimates",PETSC_DECIDE,PETSC_DECIDE,300,300,&lg);
216: NEPMonitorSet(nep,NEPMonitorLG,lg,(PetscErrorCode (*)(void**))PetscDrawLGDestroy);
217: }
218: PetscOptionsBool("-nep_monitor_lg_all","Monitor error estimates graphically","NEPMonitorSet",PETSC_FALSE,&flg,&set);
219: if (set && flg) {
220: NEPMonitorLGCreate(PetscObjectComm((PetscObject)nep),NULL,"Error estimates",PETSC_DECIDE,PETSC_DECIDE,300,300,&lg);
221: NEPMonitorSet(nep,NEPMonitorLGAll,lg,(PetscErrorCode (*)(void**))PetscDrawLGDestroy);
222: NEPSetTrackAll(nep,PETSC_TRUE);
223: }
224: /* -----------------------------------------------------------------------*/
226: PetscOptionsBoolGroupBegin("-nep_largest_magnitude","compute largest eigenvalues in magnitude","NEPSetWhichEigenpairs",&flg);
227: if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_MAGNITUDE); }
228: PetscOptionsBoolGroup("-nep_smallest_magnitude","compute smallest eigenvalues in magnitude","NEPSetWhichEigenpairs",&flg);
229: if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_MAGNITUDE); }
230: PetscOptionsBoolGroup("-nep_largest_real","compute largest real parts","NEPSetWhichEigenpairs",&flg);
231: if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_REAL); }
232: PetscOptionsBoolGroup("-nep_smallest_real","compute smallest real parts","NEPSetWhichEigenpairs",&flg);
233: if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_REAL); }
234: PetscOptionsBoolGroup("-nep_largest_imaginary","compute largest imaginary parts","NEPSetWhichEigenpairs",&flg);
235: if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_IMAGINARY); }
236: PetscOptionsBoolGroup("-nep_smallest_imaginary","compute smallest imaginary parts","NEPSetWhichEigenpairs",&flg);
237: if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_IMAGINARY); }
238: PetscOptionsBoolGroup("-nep_target_magnitude","compute nearest eigenvalues to target","NEPSetWhichEigenpairs",&flg);
239: if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE); }
240: PetscOptionsBoolGroup("-nep_target_real","compute eigenvalues with real parts close to target","NEPSetWhichEigenpairs",&flg);
241: if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_REAL); }
242: PetscOptionsBoolGroup("-nep_target_imaginary","compute eigenvalues with imaginary parts close to target","NEPSetWhichEigenpairs",&flg);
243: if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_IMAGINARY); }
244: PetscOptionsBoolGroupEnd("-nep_all","compute all eigenvalues in a region","NEPSetWhichEigenpairs",&flg);
245: if (flg) { NEPSetWhichEigenpairs(nep,NEP_ALL); }
247: PetscOptionsName("-nep_view","Print detailed information on solver used","NEPView",NULL);
248: PetscOptionsName("-nep_view_vectors","View computed eigenvectors","NEPVectorsView",NULL);
249: PetscOptionsName("-nep_view_values","View computed eigenvalues","NEPValuesView",NULL);
250: PetscOptionsName("-nep_converged_reason","Print reason for convergence, and number of iterations","NEPReasonView",NULL);
251: PetscOptionsName("-nep_error_absolute","Print absolute errors of each eigenpair","NEPErrorView",NULL);
252: PetscOptionsName("-nep_error_relative","Print relative errors of each eigenpair","NEPErrorView",NULL);
254: if (nep->ops->setfromoptions) {
255: (*nep->ops->setfromoptions)(PetscOptionsObject,nep);
256: }
257: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)nep);
258: PetscOptionsEnd();
260: if (!nep->V) { NEPGetBV(nep,&nep->V); }
261: BVSetFromOptions(nep->V);
262: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
263: RGSetFromOptions(nep->rg);
264: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
265: DSSetFromOptions(nep->ds);
266: if (!nep->refineksp) { NEPRefineGetKSP(nep,&nep->refineksp); }
267: KSPSetFromOptions(nep->refineksp);
268: return(0);
269: }
273: /*@
274: NEPGetTolerances - Gets the tolerance and maximum iteration count used
275: by the NEP convergence tests.
277: Not Collective
279: Input Parameter:
280: . nep - the nonlinear eigensolver context
282: Output Parameters:
283: + tol - the convergence tolerance
284: - maxits - maximum number of iterations
286: Notes:
287: The user can specify NULL for any parameter that is not needed.
289: Level: intermediate
291: .seealso: NEPSetTolerances()
292: @*/
293: PetscErrorCode NEPGetTolerances(NEP nep,PetscReal *tol,PetscInt *maxits)294: {
297: if (tol) *tol = nep->tol;
298: if (maxits) *maxits = nep->max_it;
299: return(0);
300: }
304: /*@
305: NEPSetTolerances - Sets the tolerance and maximum iteration count used
306: by the NEP convergence tests.
308: Logically Collective on NEP310: Input Parameters:
311: + nep - the nonlinear eigensolver context
312: . tol - the convergence tolerance
313: - maxits - maximum number of iterations to use
315: Options Database Keys:
316: + -nep_tol <tol> - Sets the convergence tolerance
317: - -nep_max_it <maxits> - Sets the maximum number of iterations allowed
319: Notes:
320: Use PETSC_DEFAULT for either argument to assign a reasonably good value.
322: Level: intermediate
324: .seealso: NEPGetTolerances()
325: @*/
326: PetscErrorCode NEPSetTolerances(NEP nep,PetscReal tol,PetscInt maxits)327: {
332: if (tol == PETSC_DEFAULT) {
333: nep->tol = PETSC_DEFAULT;
334: nep->state = NEP_STATE_INITIAL;
335: } else {
336: if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
337: nep->tol = tol;
338: }
339: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
340: nep->max_it = 0;
341: nep->state = NEP_STATE_INITIAL;
342: } else {
343: if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
344: nep->max_it = maxits;
345: }
346: return(0);
347: }
351: /*@
352: NEPGetDimensions - Gets the number of eigenvalues to compute
353: and the dimension of the subspace.
355: Not Collective
357: Input Parameter:
358: . nep - the nonlinear eigensolver context
360: Output Parameters:
361: + nev - number of eigenvalues to compute
362: . ncv - the maximum dimension of the subspace to be used by the solver
363: - mpd - the maximum dimension allowed for the projected problem
365: Notes:
366: The user can specify NULL for any parameter that is not needed.
368: Level: intermediate
370: .seealso: NEPSetDimensions()
371: @*/
372: PetscErrorCode NEPGetDimensions(NEP nep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)373: {
376: if (nev) *nev = nep->nev;
377: if (ncv) *ncv = nep->ncv;
378: if (mpd) *mpd = nep->mpd;
379: return(0);
380: }
384: /*@
385: NEPSetDimensions - Sets the number of eigenvalues to compute
386: and the dimension of the subspace.
388: Logically Collective on NEP390: Input Parameters:
391: + nep - the nonlinear eigensolver context
392: . nev - number of eigenvalues to compute
393: . ncv - the maximum dimension of the subspace to be used by the solver
394: - mpd - the maximum dimension allowed for the projected problem
396: Options Database Keys:
397: + -nep_nev <nev> - Sets the number of eigenvalues
398: . -nep_ncv <ncv> - Sets the dimension of the subspace
399: - -nep_mpd <mpd> - Sets the maximum projected dimension
401: Notes:
402: Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
403: dependent on the solution method.
405: The parameters ncv and mpd are intimately related, so that the user is advised
406: to set one of them at most. Normal usage is that
407: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
408: (b) in cases where nev is large, the user sets mpd.
410: The value of ncv should always be between nev and (nev+mpd), typically
411: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
412: a smaller value should be used.
414: Level: intermediate
416: .seealso: NEPGetDimensions()
417: @*/
418: PetscErrorCode NEPSetDimensions(NEP nep,PetscInt nev,PetscInt ncv,PetscInt mpd)419: {
425: if (nev<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
426: nep->nev = nev;
427: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
428: nep->ncv = 0;
429: } else {
430: if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
431: nep->ncv = ncv;
432: }
433: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
434: nep->mpd = 0;
435: } else {
436: if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
437: nep->mpd = mpd;
438: }
439: nep->state = NEP_STATE_INITIAL;
440: return(0);
441: }
445: /*@
446: NEPSetWhichEigenpairs - Specifies which portion of the spectrum is
447: to be sought.
449: Logically Collective on NEP451: Input Parameters:
452: + nep - eigensolver context obtained from NEPCreate()
453: - which - the portion of the spectrum to be sought
455: Possible values:
456: The parameter 'which' can have one of these values
458: + NEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
459: . NEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
460: . NEP_LARGEST_REAL - largest real parts
461: . NEP_SMALLEST_REAL - smallest real parts
462: . NEP_LARGEST_IMAGINARY - largest imaginary parts
463: . NEP_SMALLEST_IMAGINARY - smallest imaginary parts
464: . NEP_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
465: . NEP_TARGET_REAL - eigenvalues with real part closest to target
466: . NEP_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
467: . NEP_ALL - all eigenvalues contained in a given region
468: - NEP_WHICH_USER - user defined ordering set with NEPSetEigenvalueComparison()
470: Options Database Keys:
471: + -nep_largest_magnitude - Sets largest eigenvalues in magnitude
472: . -nep_smallest_magnitude - Sets smallest eigenvalues in magnitude
473: . -nep_largest_real - Sets largest real parts
474: . -nep_smallest_real - Sets smallest real parts
475: . -nep_largest_imaginary - Sets largest imaginary parts
476: . -nep_smallest_imaginary - Sets smallest imaginary parts
477: . -nep_target_magnitude - Sets eigenvalues closest to target
478: . -nep_target_real - Sets real parts closest to target
479: . -nep_target_imaginary - Sets imaginary parts closest to target
480: - -nep_all - Sets all eigenvalues in a region
482: Notes:
483: Not all eigensolvers implemented in NEP account for all the possible values
484: stated above. If SLEPc is compiled for real numbers NEP_LARGEST_IMAGINARY
485: and NEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
486: for eigenvalue selection.
488: The target is a scalar value provided with NEPSetTarget().
490: NEP_ALL is intended for use in the context of the CISS solver for
491: computing all eigenvalues in a region.
493: Level: intermediate
495: .seealso: NEPGetWhichEigenpairs(), NEPSetTarget(), NEPSetEigenvalueComparison(), NEPWhich496: @*/
497: PetscErrorCode NEPSetWhichEigenpairs(NEP nep,NEPWhich which)498: {
502: switch (which) {
503: case NEP_LARGEST_MAGNITUDE:
504: case NEP_SMALLEST_MAGNITUDE:
505: case NEP_LARGEST_REAL:
506: case NEP_SMALLEST_REAL:
507: case NEP_LARGEST_IMAGINARY:
508: case NEP_SMALLEST_IMAGINARY:
509: case NEP_TARGET_MAGNITUDE:
510: case NEP_TARGET_REAL:
511: #if defined(PETSC_USE_COMPLEX)
512: case NEP_TARGET_IMAGINARY:
513: #endif
514: case EPS_ALL:
515: case NEP_WHICH_USER:
516: if (nep->which != which) {
517: nep->state = NEP_STATE_INITIAL;
518: nep->which = which;
519: }
520: break;
521: default:522: SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
523: }
524: return(0);
525: }
529: /*@
530: NEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
531: sought.
533: Not Collective
535: Input Parameter:
536: . nep - eigensolver context obtained from NEPCreate()
538: Output Parameter:
539: . which - the portion of the spectrum to be sought
541: Notes:
542: See NEPSetWhichEigenpairs() for possible values of 'which'.
544: Level: intermediate
546: .seealso: NEPSetWhichEigenpairs(), NEPWhich547: @*/
548: PetscErrorCode NEPGetWhichEigenpairs(NEP nep,NEPWhich *which)549: {
553: *which = nep->which;
554: return(0);
555: }
559: /*@C
560: NEPSetEigenvalueComparison - Specifies the eigenvalue comparison function
561: when NEPSetWhichEigenpairs() is set to NEP_WHICH_USER.
563: Logically Collective on NEP565: Input Parameters:
566: + pep - eigensolver context obtained from NEPCreate()
567: . func - a pointer to the comparison function
568: - ctx - a context pointer (the last parameter to the comparison function)
570: Calling Sequence of func:
571: $ func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
573: + ar - real part of the 1st eigenvalue
574: . ai - imaginary part of the 1st eigenvalue
575: . br - real part of the 2nd eigenvalue
576: . bi - imaginary part of the 2nd eigenvalue
577: . res - result of comparison
578: - ctx - optional context, as set by NEPSetEigenvalueComparison()
580: Note:
581: The returning parameter 'res' can be
582: + negative - if the 1st eigenvalue is preferred to the 2st one
583: . zero - if both eigenvalues are equally preferred
584: - positive - if the 2st eigenvalue is preferred to the 1st one
586: Level: advanced
588: .seealso: NEPSetWhichEigenpairs(), NEPWhich589: @*/
590: PetscErrorCode NEPSetEigenvalueComparison(NEP pep,PetscErrorCode (*func)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)591: {
594: pep->sc->comparison = func;
595: pep->sc->comparisonctx = ctx;
596: pep->which = NEP_WHICH_USER;
597: return(0);
598: }
602: /*@C
603: NEPSetConvergenceTestFunction - Sets a function to compute the error estimate
604: used in the convergence test.
606: Logically Collective on NEP608: Input Parameters:
609: + nep - nonlinear eigensolver context obtained from NEPCreate()
610: . func - a pointer to the convergence test function
611: . ctx - context for private data for the convergence routine (may be null)
612: - destroy - a routine for destroying the context (may be null)
614: Calling Sequence of func:
615: $ func(NEP nep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
617: + nep - nonlinear eigensolver context obtained from NEPCreate()
618: . eigr - real part of the eigenvalue
619: . eigi - imaginary part of the eigenvalue
620: . res - residual norm associated to the eigenpair
621: . errest - (output) computed error estimate
622: - ctx - optional context, as set by NEPSetConvergenceTestFunction()
624: Note:
625: If the error estimate returned by the convergence test function is less than
626: the tolerance, then the eigenvalue is accepted as converged.
628: Level: advanced
630: .seealso: NEPSetConvergenceTest(), NEPSetTolerances()
631: @*/
632: PetscErrorCode NEPSetConvergenceTestFunction(NEP nep,PetscErrorCode (*func)(NEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))633: {
638: if (nep->convergeddestroy) {
639: (*nep->convergeddestroy)(nep->convergedctx);
640: }
641: nep->converged = func;
642: nep->convergeddestroy = destroy;
643: nep->convergedctx = ctx;
644: if (func == NEPConvergedRelative) nep->conv = NEP_CONV_REL;
645: else if (func == NEPConvergedNorm) nep->conv = NEP_CONV_NORM;
646: else if (func == NEPConvergedAbsolute) nep->conv = NEP_CONV_ABS;
647: else nep->conv = NEP_CONV_USER;
648: return(0);
649: }
653: /*@
654: NEPSetConvergenceTest - Specifies how to compute the error estimate
655: used in the convergence test.
657: Logically Collective on NEP659: Input Parameters:
660: + nep - nonlinear eigensolver context obtained from NEPCreate()
661: - conv - the type of convergence test
663: Options Database Keys:
664: + -nep_conv_abs - Sets the absolute convergence test
665: . -nep_conv_rel - Sets the convergence test relative to the eigenvalue
666: - -nep_conv_user - Selects the user-defined convergence test
668: Note:
669: The parameter 'conv' can have one of these values
670: + NEP_CONV_ABS - absolute error ||r||
671: . NEP_CONV_REL - error relative to the eigenvalue l, ||r||/|l|
672: . NEP_CONV_NORM - error relative matrix norms, ||r||/sum_i(|f_i(l)|*||A_i||)
673: - NEP_CONV_USER - function set by NEPSetConvergenceTestFunction()
675: Level: intermediate
677: .seealso: NEPGetConvergenceTest(), NEPSetConvergenceTestFunction(), NEPSetStoppingTest(), NEPConv678: @*/
679: PetscErrorCode NEPSetConvergenceTest(NEP nep,NEPConv conv)680: {
684: switch (conv) {
685: case NEP_CONV_ABS: nep->converged = NEPConvergedAbsolute; break;
686: case NEP_CONV_REL: nep->converged = NEPConvergedRelative; break;
687: case NEP_CONV_NORM: nep->converged = NEPConvergedNorm; break;
688: case NEP_CONV_USER: break;
689: default:690: SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
691: }
692: nep->conv = conv;
693: return(0);
694: }
698: /*@
699: NEPGetConvergenceTest - Gets the method used to compute the error estimate
700: used in the convergence test.
702: Not Collective
704: Input Parameters:
705: . nep - nonlinear eigensolver context obtained from NEPCreate()
707: Output Parameters:
708: . conv - the type of convergence test
710: Level: intermediate
712: .seealso: NEPSetConvergenceTest(), NEPConv713: @*/
714: PetscErrorCode NEPGetConvergenceTest(NEP nep,NEPConv *conv)715: {
719: *conv = nep->conv;
720: return(0);
721: }
725: /*@C
726: NEPSetStoppingTestFunction - Sets a function to decide when to stop the outer
727: iteration of the eigensolver.
729: Logically Collective on NEP731: Input Parameters:
732: + nep - nonlinear eigensolver context obtained from NEPCreate()
733: . func - pointer to the stopping test function
734: . ctx - context for private data for the stopping routine (may be null)
735: - destroy - a routine for destroying the context (may be null)
737: Calling Sequence of func:
738: $ func(NEP nep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,NEPConvergedReason *reason,void *ctx)
740: + nep - nonlinear eigensolver context obtained from NEPCreate()
741: . its - current number of iterations
742: . max_it - maximum number of iterations
743: . nconv - number of currently converged eigenpairs
744: . nev - number of requested eigenpairs
745: . reason - (output) result of the stopping test
746: - ctx - optional context, as set by NEPSetStoppingTestFunction()
748: Note:
749: Normal usage is to first call the default routine NEPStoppingBasic() and then
750: set reason to NEP_CONVERGED_USER if some user-defined conditions have been
751: met. To let the eigensolver continue iterating, the result must be left as
752: NEP_CONVERGED_ITERATING.
754: Level: advanced
756: .seealso: NEPSetStoppingTest(), NEPStoppingBasic()
757: @*/
758: PetscErrorCode NEPSetStoppingTestFunction(NEP nep,PetscErrorCode (*func)(NEP,PetscInt,PetscInt,PetscInt,PetscInt,NEPConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))759: {
764: if (nep->stoppingdestroy) {
765: (*nep->stoppingdestroy)(nep->stoppingctx);
766: }
767: nep->stopping = func;
768: nep->stoppingdestroy = destroy;
769: nep->stoppingctx = ctx;
770: if (func == NEPStoppingBasic) nep->stop = NEP_STOP_BASIC;
771: else nep->stop = NEP_STOP_USER;
772: return(0);
773: }
777: /*@
778: NEPSetStoppingTest - Specifies how to decide the termination of the outer
779: loop of the eigensolver.
781: Logically Collective on NEP783: Input Parameters:
784: + nep - nonlinear eigensolver context obtained from NEPCreate()
785: - stop - the type of stopping test
787: Options Database Keys:
788: + -nep_stop_basic - Sets the default stopping test
789: - -nep_stop_user - Selects the user-defined stopping test
791: Note:
792: The parameter 'stop' can have one of these values
793: + NEP_STOP_BASIC - default stopping test
794: - NEP_STOP_USER - function set by NEPSetStoppingTestFunction()
796: Level: advanced
798: .seealso: NEPGetStoppingTest(), NEPSetStoppingTestFunction(), NEPSetConvergenceTest(), NEPStop799: @*/
800: PetscErrorCode NEPSetStoppingTest(NEP nep,NEPStop stop)801: {
805: switch (stop) {
806: case NEP_STOP_BASIC: nep->stopping = NEPStoppingBasic; break;
807: case NEP_STOP_USER: break;
808: default:809: SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
810: }
811: nep->stop = stop;
812: return(0);
813: }
817: /*@
818: NEPGetStoppingTest - Gets the method used to decide the termination of the outer
819: loop of the eigensolver.
821: Not Collective
823: Input Parameters:
824: . nep - nonlinear eigensolver context obtained from NEPCreate()
826: Output Parameters:
827: . stop - the type of stopping test
829: Level: advanced
831: .seealso: NEPSetStoppingTest(), NEPStop832: @*/
833: PetscErrorCode NEPGetStoppingTest(NEP nep,NEPStop *stop)834: {
838: *stop = nep->stop;
839: return(0);
840: }
844: /*@
845: NEPSetTrackAll - Specifies if the solver must compute the residual of all
846: approximate eigenpairs or not.
848: Logically Collective on NEP850: Input Parameters:
851: + nep - the eigensolver context
852: - trackall - whether compute all residuals or not
854: Notes:
855: If the user sets trackall=PETSC_TRUE then the solver explicitly computes
856: the residual for each eigenpair approximation. Computing the residual is
857: usually an expensive operation and solvers commonly compute the associated
858: residual to the first unconverged eigenpair.
859: The options '-nep_monitor_all' and '-nep_monitor_lg_all' automatically
860: activate this option.
862: Level: developer
864: .seealso: NEPGetTrackAll()
865: @*/
866: PetscErrorCode NEPSetTrackAll(NEP nep,PetscBool trackall)867: {
871: nep->trackall = trackall;
872: return(0);
873: }
877: /*@
878: NEPGetTrackAll - Returns the flag indicating whether all residual norms must
879: be computed or not.
881: Not Collective
883: Input Parameter:
884: . nep - the eigensolver context
886: Output Parameter:
887: . trackall - the returned flag
889: Level: developer
891: .seealso: NEPSetTrackAll()
892: @*/
893: PetscErrorCode NEPGetTrackAll(NEP nep,PetscBool *trackall)894: {
898: *trackall = nep->trackall;
899: return(0);
900: }
904: /*@
905: NEPSetRefine - Specifies the refinement type (and options) to be used
906: after the solve.
908: Logically Collective on NEP910: Input Parameters:
911: + nep - the nonlinear eigensolver context
912: . refine - refinement type
913: . npart - number of partitions of the communicator
914: . tol - the convergence tolerance
915: . its - maximum number of refinement iterations
916: - scheme - which scheme to be used for solving the involved linear systems
918: Options Database Keys:
919: + -nep_refine <type> - refinement type, one of <none,simple,multiple>
920: . -nep_refine_partitions <n> - the number of partitions
921: . -nep_refine_tol <tol> - the tolerance
922: . -nep_refine_its <its> - number of iterations
923: - -nep_refine_scheme - to set the scheme for the linear solves
925: Notes:
926: By default, iterative refinement is disabled, since it may be very
927: costly. There are two possible refinement strategies: simple and multiple.
928: The simple approach performs iterative refinement on each of the
929: converged eigenpairs individually, whereas the multiple strategy works
930: with the invariant pair as a whole, refining all eigenpairs simultaneously.
931: The latter may be required for the case of multiple eigenvalues.
933: In some cases, especially when using direct solvers within the
934: iterative refinement method, it may be helpful for improved scalability
935: to split the communicator in several partitions. The npart parameter
936: indicates how many partitions to use (defaults to 1).
938: The tol and its parameters specify the stopping criterion. In the simple
939: method, refinement continues until the residual of each eigenpair is
940: below the tolerance (tol defaults to the NEP tol, but may be set to a
941: different value). In contrast, the multiple method simply performs its
942: refinement iterations (just one by default).
944: The scheme argument is used to change the way in which linear systems are
945: solved. Possible choices are: explicit, mixed block elimination (MBE),
946: and Schur complement.
948: Level: intermediate
950: .seealso: NEPGetRefine()
951: @*/
952: PetscErrorCode NEPSetRefine(NEP nep,NEPRefine refine,PetscInt npart,PetscReal tol,PetscInt its,NEPRefineScheme scheme)953: {
955: PetscMPIInt size;
964: nep->refine = refine;
965: if (refine) { /* process parameters only if not REFINE_NONE */
966: if (npart == PETSC_DEFAULT || npart == PETSC_DECIDE) {
967: nep->npart = 1;
968: } else {
969: MPI_Comm_size(PetscObjectComm((PetscObject)nep),&size);
970: if (npart<1 || npart>size) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of npart");
971: nep->npart = npart;
972: }
973: if (tol == PETSC_DEFAULT || tol == PETSC_DECIDE) {
974: nep->rtol = PetscMax(nep->tol/1000,PETSC_MACHINE_EPSILON);
975: } else {
976: if (tol<=0.0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
977: nep->rtol = tol;
978: }
979: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) {
980: nep->rits = PETSC_DEFAULT;
981: } else {
982: if (its<0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be >= 0");
983: nep->rits = its;
984: }
985: nep->scheme = scheme;
986: }
987: nep->state = NEP_STATE_INITIAL;
988: return(0);
989: }
993: /*@
994: NEPGetRefine - Gets the refinement strategy used by the NEP object, and the
995: associated parameters.
997: Not Collective
999: Input Parameter:
1000: . nep - the nonlinear eigensolver context
1002: Output Parameters:
1003: + refine - refinement type
1004: . npart - number of partitions of the communicator
1005: . tol - the convergence tolerance
1006: - its - maximum number of refinement iterations
1007: - scheme - the scheme used for solving linear systems
1009: Level: intermediate
1011: Note:
1012: The user can specify NULL for any parameter that is not needed.
1014: .seealso: NEPSetRefine()
1015: @*/
1016: PetscErrorCode NEPGetRefine(NEP nep,NEPRefine *refine,PetscInt *npart,PetscReal *tol,PetscInt *its,NEPRefineScheme *scheme)1017: {
1020: if (refine) *refine = nep->refine;
1021: if (npart) *npart = nep->npart;
1022: if (tol) *tol = nep->rtol;
1023: if (its) *its = nep->rits;
1024: if (scheme) *scheme = nep->scheme;
1025: return(0);
1026: }
1030: /*@C
1031: NEPSetOptionsPrefix - Sets the prefix used for searching for all
1032: NEP options in the database.
1034: Logically Collective on NEP1036: Input Parameters:
1037: + nep - the nonlinear eigensolver context
1038: - prefix - the prefix string to prepend to all NEP option requests
1040: Notes:
1041: A hyphen (-) must NOT be given at the beginning of the prefix name.
1042: The first character of all runtime options is AUTOMATICALLY the
1043: hyphen.
1045: For example, to distinguish between the runtime options for two
1046: different NEP contexts, one could call
1047: .vb
1048: NEPSetOptionsPrefix(nep1,"neig1_")
1049: NEPSetOptionsPrefix(nep2,"neig2_")
1050: .ve
1052: Level: advanced
1054: .seealso: NEPAppendOptionsPrefix(), NEPGetOptionsPrefix()
1055: @*/
1056: PetscErrorCode NEPSetOptionsPrefix(NEP nep,const char *prefix)1057: {
1062: if (!nep->V) { NEPGetBV(nep,&nep->V); }
1063: BVSetOptionsPrefix(nep->V,prefix);
1064: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
1065: DSSetOptionsPrefix(nep->ds,prefix);
1066: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
1067: RGSetOptionsPrefix(nep->rg,prefix);
1068: PetscObjectSetOptionsPrefix((PetscObject)nep,prefix);
1069: return(0);
1070: }
1074: /*@C
1075: NEPAppendOptionsPrefix - Appends to the prefix used for searching for all
1076: NEP options in the database.
1078: Logically Collective on NEP1080: Input Parameters:
1081: + nep - the nonlinear eigensolver context
1082: - prefix - the prefix string to prepend to all NEP option requests
1084: Notes:
1085: A hyphen (-) must NOT be given at the beginning of the prefix name.
1086: The first character of all runtime options is AUTOMATICALLY the hyphen.
1088: Level: advanced
1090: .seealso: NEPSetOptionsPrefix(), NEPGetOptionsPrefix()
1091: @*/
1092: PetscErrorCode NEPAppendOptionsPrefix(NEP nep,const char *prefix)1093: {
1098: if (!nep->V) { NEPGetBV(nep,&nep->V); }
1099: BVSetOptionsPrefix(nep->V,prefix);
1100: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
1101: DSSetOptionsPrefix(nep->ds,prefix);
1102: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
1103: RGSetOptionsPrefix(nep->rg,prefix);
1104: PetscObjectAppendOptionsPrefix((PetscObject)nep,prefix);
1105: return(0);
1106: }
1110: /*@C
1111: NEPGetOptionsPrefix - Gets the prefix used for searching for all
1112: NEP options in the database.
1114: Not Collective
1116: Input Parameters:
1117: . nep - the nonlinear eigensolver context
1119: Output Parameters:
1120: . prefix - pointer to the prefix string used is returned
1122: Note:
1123: On the Fortran side, the user should pass in a string 'prefix' of
1124: sufficient length to hold the prefix.
1126: Level: advanced
1128: .seealso: NEPSetOptionsPrefix(), NEPAppendOptionsPrefix()
1129: @*/
1130: PetscErrorCode NEPGetOptionsPrefix(NEP nep,const char *prefix[])1131: {
1137: PetscObjectGetOptionsPrefix((PetscObject)nep,prefix);
1138: return(0);
1139: }