Actual source code: nepopts.c

slepc-3.5.2 2014-10-10
Report Typos and Errors
  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-2014, 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,flg1,flg2,flg3,flg4,flg5;
 49:   PetscReal        r1,r2,r3;
 50:   PetscScalar      s;
 51:   PetscInt         i,j,k;
 52:   PetscViewer      monviewer;
 53:   SlepcConvMonitor ctx;

 57:   if (!NEPRegisterAllCalled) { NEPRegisterAll(); }
 58:   PetscObjectOptionsBegin((PetscObject)nep);
 59:     PetscOptionsFList("-nep_type","Nonlinear Eigenvalue Problem method","NEPSetType",NEPList,(char*)(((PetscObject)nep)->type_name?((PetscObject)nep)->type_name:NEPRII),type,256,&flg);
 60:     if (flg) {
 61:       NEPSetType(nep,type);
 62:     } else if (!((PetscObject)nep)->type_name) {
 63:       NEPSetType(nep,NEPRII);
 64:     }

 66:     PetscOptionsEnum("-nep_refine","Iterative refinement method","NEPSetRefine",NEPRefineTypes,(PetscEnum)nep->refine,(PetscEnum*)&nep->refine,NULL);

 68:     r1 = nep->reftol;
 69:     PetscOptionsReal("-nep_refine_tol","Tolerance for iterative refinement","NEPSetRefine",nep->reftol,&r1,&flg1);
 70:     j = nep->rits;
 71:     PetscOptionsInt("-nep_refine_its","Maximum number of iterations for iterative refinement","NEPSetRefine",nep->rits,&j,&flg2);
 72:     if (flg1 || flg2) {
 73:       NEPSetRefine(nep,nep->refine,r1,j);
 74:     }

 76:     i = nep->max_it? nep->max_it: PETSC_DEFAULT;
 77:     PetscOptionsInt("-nep_max_it","Maximum number of iterations","NEPSetTolerances",nep->max_it,&i,&flg1);
 78:     j = nep->max_funcs? nep->max_funcs: PETSC_DEFAULT;
 79:     PetscOptionsInt("-nep_max_funcs","Maximum number of function evaluations","NEPSetTolerances",nep->max_funcs,&j,&flg2);
 80:     r1 = nep->abstol;
 81:     PetscOptionsReal("-nep_atol","Absolute tolerance for residual norm","NEPSetTolerances",nep->abstol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:nep->abstol,&r1,&flg3);
 82:     r2 = nep->rtol;
 83:     PetscOptionsReal("-nep_rtol","Relative tolerance for residual norm","NEPSetTolerances",nep->rtol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:nep->rtol,&r2,&flg4);
 84:     r3 = nep->stol;
 85:     PetscOptionsReal("-nep_stol","Relative tolerance for step length","NEPSetTolerances",nep->stol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:nep->stol,&r3,&flg5);
 86:     if (flg1 || flg2 || flg3 || flg4 || flg5) {
 87:       NEPSetTolerances(nep,r1,r2,r3,i,j);
 88:     }

 90:     flg  = PETSC_FALSE;
 91:     PetscOptionsBool("-nep_convergence_default","Default (relative error) convergence test","NEPSetConvergenceTest",flg,&flg,NULL);
 92:     if (flg) {
 93:       NEPSetConvergenceTest(nep,NEPConvergedDefault,NULL,NULL);
 94:     }

 96:     i = nep->nev;
 97:     PetscOptionsInt("-nep_nev","Number of eigenvalues to compute","NEPSetDimensions",nep->nev,&i,&flg1);
 98:     j = nep->ncv? nep->ncv: PETSC_DEFAULT;
 99:     PetscOptionsInt("-nep_ncv","Number of basis vectors","NEPSetDimensions",nep->ncv,&j,&flg2);
100:     k = nep->mpd? nep->mpd: PETSC_DEFAULT;
101:     PetscOptionsInt("-nep_mpd","Maximum dimension of projected problem","NEPSetDimensions",nep->mpd,&k,&flg3);
102:     if (flg1 || flg2 || flg3) {
103:       NEPSetDimensions(nep,i,j,k);
104:     }

106:     i = 0;
107:     PetscOptionsInt("-nep_lag_preconditioner","Interval to rebuild preconditioner","NEPSetLagPreconditioner",nep->lag,&i,&flg);
108:     if (flg) { NEPSetLagPreconditioner(nep,i); }

110:     PetscOptionsBool("-nep_const_correction_tol","Constant correction tolerance for the linear solver","NEPSetConstCorrectionTol",nep->cctol,&nep->cctol,NULL);

112:     PetscOptionsScalar("-nep_target","Value of the target","NEPSetTarget",nep->target,&s,&flg);
113:     if (flg) {
114:       NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE);
115:       NEPSetTarget(nep,s);
116:     }

118:     /* -----------------------------------------------------------------------*/
119:     /*
120:       Cancels all monitors hardwired into code before call to NEPSetFromOptions()
121:     */
122:     flg = PETSC_FALSE;
123:     PetscOptionsBool("-nep_monitor_cancel","Remove any hardwired monitor routines","NEPMonitorCancel",flg,&flg,NULL);
124:     if (flg) {
125:       NEPMonitorCancel(nep);
126:     }
127:     /*
128:       Prints approximate eigenvalues and error estimates at each iteration
129:     */
130:     PetscOptionsString("-nep_monitor","Monitor first unconverged approximate eigenvalue and error estimate","NEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
131:     if (flg) {
132:       PetscViewerASCIIOpen(PetscObjectComm((PetscObject)nep),monfilename,&monviewer);
133:       NEPMonitorSet(nep,NEPMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
134:     }
135:     PetscOptionsString("-nep_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","NEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
136:     if (flg) {
137:       PetscNew(&ctx);
138:       PetscViewerASCIIOpen(PetscObjectComm((PetscObject)nep),monfilename,&ctx->viewer);
139:       NEPMonitorSet(nep,NEPMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
140:     }
141:     PetscOptionsString("-nep_monitor_all","Monitor approximate eigenvalues and error estimates","NEPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
142:     if (flg) {
143:       PetscViewerASCIIOpen(PetscObjectComm((PetscObject)nep),monfilename,&monviewer);
144:       NEPMonitorSet(nep,NEPMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
145:       NEPSetTrackAll(nep,PETSC_TRUE);
146:     }
147:     flg = PETSC_FALSE;
148:     PetscOptionsBool("-nep_monitor_lg","Monitor first unconverged approximate error estimate graphically","NEPMonitorSet",flg,&flg,NULL);
149:     if (flg) {
150:       NEPMonitorSet(nep,NEPMonitorLG,NULL,NULL);
151:     }
152:     flg = PETSC_FALSE;
153:     PetscOptionsBool("-nep_monitor_lg_all","Monitor error estimates graphically","NEPMonitorSet",flg,&flg,NULL);
154:     if (flg) {
155:       NEPMonitorSet(nep,NEPMonitorLGAll,NULL,NULL);
156:       NEPSetTrackAll(nep,PETSC_TRUE);
157:     }
158:   /* -----------------------------------------------------------------------*/

160:     PetscOptionsBoolGroupBegin("-nep_largest_magnitude","compute largest eigenvalues in magnitude","NEPSetWhichEigenpairs",&flg);
161:     if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_MAGNITUDE); }
162:     PetscOptionsBoolGroup("-nep_smallest_magnitude","compute smallest eigenvalues in magnitude","NEPSetWhichEigenpairs",&flg);
163:     if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_MAGNITUDE); }
164:     PetscOptionsBoolGroup("-nep_largest_real","compute largest real parts","NEPSetWhichEigenpairs",&flg);
165:     if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_REAL); }
166:     PetscOptionsBoolGroup("-nep_smallest_real","compute smallest real parts","NEPSetWhichEigenpairs",&flg);
167:     if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_REAL); }
168:     PetscOptionsBoolGroup("-nep_largest_imaginary","compute largest imaginary parts","NEPSetWhichEigenpairs",&flg);
169:     if (flg) { NEPSetWhichEigenpairs(nep,NEP_LARGEST_IMAGINARY); }
170:     PetscOptionsBoolGroup("-nep_smallest_imaginary","compute smallest imaginary parts","NEPSetWhichEigenpairs",&flg);
171:     if (flg) { NEPSetWhichEigenpairs(nep,NEP_SMALLEST_IMAGINARY); }
172:     PetscOptionsBoolGroup("-nep_target_magnitude","compute nearest eigenvalues to target","NEPSetWhichEigenpairs",&flg);
173:     if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE); }
174:     PetscOptionsBoolGroup("-nep_target_real","compute eigenvalues with real parts close to target","NEPSetWhichEigenpairs",&flg);
175:     if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_REAL); }
176:     PetscOptionsBoolGroupEnd("-nep_target_imaginary","compute eigenvalues with imaginary parts close to target","NEPSetWhichEigenpairs",&flg);
177:     if (flg) { NEPSetWhichEigenpairs(nep,NEP_TARGET_IMAGINARY); }

179:     PetscOptionsName("-nep_view","Print detailed information on solver used","NEPView",0);
180:     PetscOptionsName("-nep_plot_eigs","Make a plot of the computed eigenvalues","NEPSolve",0);

182:     if (nep->ops->setfromoptions) {
183:       (*nep->ops->setfromoptions)(nep);
184:     }
185:     PetscObjectProcessOptionsHandlers((PetscObject)nep);
186:   PetscOptionsEnd();

188:   if (!nep->V) { NEPGetBV(nep,&nep->V); }
189:   BVSetFromOptions(nep->V);
190:   if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
191:   RGSetFromOptions(nep->rg);
192:   if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
193:   DSSetFromOptions(nep->ds);
194:   if (!nep->ksp) { NEPGetKSP(nep,&nep->ksp); }
195:   KSPSetOperators(nep->ksp,nep->function,nep->function_pre);
196:   KSPSetFromOptions(nep->ksp);
197:   PetscRandomSetFromOptions(nep->rand);
198:   return(0);
199: }

203: /*@
204:    NEPGetTolerances - Gets the tolerance and maximum iteration count used
205:    by the NEP convergence tests.

207:    Not Collective

209:    Input Parameter:
210: .  nep - the nonlinear eigensolver context

212:    Output Parameters:
213: +  abstol - absolute convergence tolerance
214: .  rtol   - relative convergence tolerance
215: .  stol   - convergence tolerance in terms of the norm of the change in the
216:            solution between steps, || delta x || < stol*|| x ||
217: .  maxit  - maximum number of iterations
218: -  maxf   - maximum number of function evaluations

220:    Notes:
221:    The user can specify NULL for any parameter that is not needed.

223:    Level: intermediate

225: .seealso: NEPSetTolerances()
226: @*/
227: PetscErrorCode NEPGetTolerances(NEP nep,PetscReal *abstol,PetscReal *rtol,PetscReal *stol,PetscInt *maxit,PetscInt *maxf)
228: {
231:   if (abstol) *abstol = nep->abstol;
232:   if (rtol)   *rtol   = nep->rtol;
233:   if (stol)   *stol   = nep->stol;
234:   if (maxit)  *maxit  = nep->max_it;
235:   if (maxf)   *maxf   = nep->max_funcs;
236:   return(0);
237: }

241: /*@
242:    NEPSetTolerances - Sets various parameters used in convergence tests.

244:    Logically Collective on NEP

246:    Input Parameters:
247: +  nep    - the nonlinear eigensolver context
248: .  abstol - absolute convergence tolerance
249: .  rtol   - relative convergence tolerance
250: .  stol   - convergence tolerance in terms of the norm of the change in the
251:             solution between steps, || delta x || < stol*|| x ||
252: .  maxit  - maximum number of iterations
253: -  maxf   - maximum number of function evaluations

255:    Options Database Keys:
256: +    -nep_atol <abstol> - Sets abstol
257: .    -nep_rtol <rtol> - Sets rtol
258: .    -nep_stol <stol> - Sets stol
259: .    -nep_max_it <maxit> - Sets maxit
260: -    -nep_max_funcs <maxf> - Sets maxf

262:    Notes:
263:    Use PETSC_DEFAULT for either argument to assign a reasonably good value.

265:    Level: intermediate

267: .seealso: NEPGetTolerances()
268: @*/
269: PetscErrorCode NEPSetTolerances(NEP nep,PetscReal abstol,PetscReal rtol,PetscReal stol,PetscInt maxit,PetscInt maxf)
270: {
278:   if (abstol == PETSC_DEFAULT) {
279:     nep->abstol = PETSC_DEFAULT;
280:   } else {
281:     if (abstol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %g must be non-negative",(double)abstol);
282:     nep->abstol = abstol;
283:   }
284:   if (rtol == PETSC_DEFAULT) {
285:     nep->rtol = PETSC_DEFAULT;
286:   } else {
287:     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",(double)rtol);
288:     nep->rtol = rtol;
289:   }
290:   if (stol == PETSC_DEFAULT) {
291:     nep->stol = PETSC_DEFAULT;
292:   } else {
293:     if (stol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Step tolerance %g must be non-negative",(double)stol);
294:     nep->stol = stol;
295:   }
296:   if (maxit == PETSC_DEFAULT || maxit == PETSC_DECIDE) {
297:     nep->max_it = 0;
298:     nep->state = NEP_STATE_INITIAL;
299:   } else {
300:     if (maxit < 0) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",maxit);
301:     nep->max_it = maxit;
302:   }
303:   if (maxf == PETSC_DEFAULT || maxf == PETSC_DECIDE) {
304:     nep->max_it = 0;
305:     nep->state = NEP_STATE_INITIAL;
306:   } else {
307:     if (maxf < 0) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of function evaluations %D must be non-negative",maxf);
308:     nep->max_funcs = maxf;
309:   }
310:   return(0);
311: }

315: /*@
316:    NEPGetDimensions - Gets the number of eigenvalues to compute
317:    and the dimension of the subspace.

319:    Not Collective

321:    Input Parameter:
322: .  nep - the nonlinear eigensolver context

324:    Output Parameters:
325: +  nev - number of eigenvalues to compute
326: .  ncv - the maximum dimension of the subspace to be used by the solver
327: -  mpd - the maximum dimension allowed for the projected problem

329:    Notes:
330:    The user can specify NULL for any parameter that is not needed.

332:    Level: intermediate

334: .seealso: NEPSetDimensions()
335: @*/
336: PetscErrorCode NEPGetDimensions(NEP nep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
337: {
340:   if (nev) *nev = nep->nev;
341:   if (ncv) *ncv = nep->ncv;
342:   if (mpd) *mpd = nep->mpd;
343:   return(0);
344: }

348: /*@
349:    NEPSetDimensions - Sets the number of eigenvalues to compute
350:    and the dimension of the subspace.

352:    Logically Collective on NEP

354:    Input Parameters:
355: +  nep - the nonlinear eigensolver context
356: .  nev - number of eigenvalues to compute
357: .  ncv - the maximum dimension of the subspace to be used by the solver
358: -  mpd - the maximum dimension allowed for the projected problem

360:    Options Database Keys:
361: +  -nep_nev <nev> - Sets the number of eigenvalues
362: .  -nep_ncv <ncv> - Sets the dimension of the subspace
363: -  -nep_mpd <mpd> - Sets the maximum projected dimension

365:    Notes:
366:    Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
367:    dependent on the solution method.

369:    The parameters ncv and mpd are intimately related, so that the user is advised
370:    to set one of them at most. Normal usage is that
371:    (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
372:    (b) in cases where nev is large, the user sets mpd.

374:    The value of ncv should always be between nev and (nev+mpd), typically
375:    ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
376:    a smaller value should be used.

378:    Level: intermediate

380: .seealso: NEPGetDimensions()
381: @*/
382: PetscErrorCode NEPSetDimensions(NEP nep,PetscInt nev,PetscInt ncv,PetscInt mpd)
383: {
389:   if (nev<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
390:   nep->nev = nev;
391:   if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
392:     nep->ncv = 0;
393:   } else {
394:     if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
395:     nep->ncv = ncv;
396:   }
397:   if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
398:     nep->mpd = 0;
399:   } else {
400:     if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
401:     nep->mpd = mpd;
402:   }
403:   nep->state = NEP_STATE_INITIAL;
404:   return(0);
405: }

409: /*@
410:     NEPSetWhichEigenpairs - Specifies which portion of the spectrum is
411:     to be sought.

413:     Logically Collective on NEP

415:     Input Parameters:
416: +   nep   - eigensolver context obtained from NEPCreate()
417: -   which - the portion of the spectrum to be sought

419:     Possible values:
420:     The parameter 'which' can have one of these values

422: +     NEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
423: .     NEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
424: .     NEP_LARGEST_REAL - largest real parts
425: .     NEP_SMALLEST_REAL - smallest real parts
426: .     NEP_LARGEST_IMAGINARY - largest imaginary parts
427: .     NEP_SMALLEST_IMAGINARY - smallest imaginary parts
428: .     NEP_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
429: .     NEP_TARGET_REAL - eigenvalues with real part closest to target
430: -     NEP_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target

432:     Options Database Keys:
433: +   -nep_largest_magnitude - Sets largest eigenvalues in magnitude
434: .   -nep_smallest_magnitude - Sets smallest eigenvalues in magnitude
435: .   -nep_largest_real - Sets largest real parts
436: .   -nep_smallest_real - Sets smallest real parts
437: .   -nep_largest_imaginary - Sets largest imaginary parts
438: .   -nep_smallest_imaginary - Sets smallest imaginary parts
439: .   -nep_target_magnitude - Sets eigenvalues closest to target
440: .   -nep_target_real - Sets real parts closest to target
441: -   -nep_target_imaginary - Sets imaginary parts closest to target

443:     Notes:
444:     Not all eigensolvers implemented in NEP account for all the possible values
445:     stated above. If SLEPc is compiled for real numbers NEP_LARGEST_IMAGINARY
446:     and NEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
447:     for eigenvalue selection.

449:     Level: intermediate

451: .seealso: NEPGetWhichEigenpairs(), NEPWhich
452: @*/
453: PetscErrorCode NEPSetWhichEigenpairs(NEP nep,NEPWhich which)
454: {
458:   switch (which) {
459:     case NEP_LARGEST_MAGNITUDE:
460:     case NEP_SMALLEST_MAGNITUDE:
461:     case NEP_LARGEST_REAL:
462:     case NEP_SMALLEST_REAL:
463:     case NEP_LARGEST_IMAGINARY:
464:     case NEP_SMALLEST_IMAGINARY:
465:     case NEP_TARGET_MAGNITUDE:
466:     case NEP_TARGET_REAL:
467: #if defined(PETSC_USE_COMPLEX)
468:     case NEP_TARGET_IMAGINARY:
469: #endif
470:       if (nep->which != which) {
471:         nep->state = NEP_STATE_INITIAL;
472:         nep->which = which;
473:       }
474:       break;
475:     default:
476:       SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
477:   }
478:   return(0);
479: }

483: /*@C
484:     NEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
485:     sought.

487:     Not Collective

489:     Input Parameter:
490: .   nep - eigensolver context obtained from NEPCreate()

492:     Output Parameter:
493: .   which - the portion of the spectrum to be sought

495:     Notes:
496:     See NEPSetWhichEigenpairs() for possible values of 'which'.

498:     Level: intermediate

500: .seealso: NEPSetWhichEigenpairs(), NEPWhich
501: @*/
502: PetscErrorCode NEPGetWhichEigenpairs(NEP nep,NEPWhich *which)
503: {
507:   *which = nep->which;
508:   return(0);
509: }

513: /*@
514:     NEPSetLagPreconditioner - Determines when the preconditioner is rebuilt in the
515:     nonlinear solve.

517:     Logically Collective on NEP

519:     Input Parameters:
520: +   nep - the NEP context
521: -   lag - 0 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is
522:           computed within the nonlinear iteration, 2 means every second time
523:           the Jacobian is built, etc.

525:     Options Database Keys:
526: .   -nep_lag_preconditioner <lag>

528:     Notes:
529:     The default is 1.
530:     The preconditioner is ALWAYS built in the first iteration of a nonlinear solve.

532:     Level: intermediate

534: .seealso: NEPGetLagPreconditioner()
535: @*/
536: PetscErrorCode NEPSetLagPreconditioner(NEP nep,PetscInt lag)
537: {
541:   if (lag<0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag must be non-negative");
542:   nep->lag = lag;
543:   return(0);
544: }

548: /*@
549:     NEPGetLagPreconditioner - Indicates how often the preconditioner is rebuilt.

551:     Not Collective

553:     Input Parameter:
554: .   nep - the NEP context

556:     Output Parameter:
557: .   lag - the lag parameter

559:     Level: intermediate

561: .seealso: NEPSetLagPreconditioner()
562: @*/
563: PetscErrorCode NEPGetLagPreconditioner(NEP nep,PetscInt *lag)
564: {
568:   *lag = nep->lag;
569:   return(0);
570: }

574: /*@
575:     NEPSetConstCorrectionTol - Sets a flag to keep the tolerance used
576:     in the linear solver constant.

578:     Logically Collective on NEP

580:     Input Parameters:
581: +   nep - the NEP context
582: -   cct - a boolean value

584:     Options Database Keys:
585: .   -nep_const_correction_tol <cct>

587:     Notes:
588:     By default, an exponentially decreasing tolerance is set in the KSP used
589:     within the nonlinear iteration, so that each Newton iteration requests
590:     better accuracy than the previous one. The constant correction tolerance
591:     flag stops this behaviour.

593:     Level: intermediate

595: .seealso: NEPGetConstCorrectionTol()
596: @*/
597: PetscErrorCode NEPSetConstCorrectionTol(NEP nep,PetscBool cct)
598: {
602:   nep->cctol = cct;
603:   return(0);
604: }

608: /*@
609:     NEPGetConstCorrectionTol - Returns the constant tolerance flag.

611:     Not Collective

613:     Input Parameter:
614: .   nep - the NEP context

616:     Output Parameter:
617: .   cct - the value of the constant tolerance flag

619:     Level: intermediate

621: .seealso: NEPSetConstCorrectionTol()
622: @*/
623: PetscErrorCode NEPGetConstCorrectionTol(NEP nep,PetscBool *cct)
624: {
628:   *cct = nep->cctol;
629:   return(0);
630: }

634: /*@C
635:     NEPSetConvergenceTest - Sets the function to be used to test convergence
636:     of the nonlinear iterative solution.

638:     Logically Collective on NEP

640:     Input Parameters:
641: +   nep     - the NEP context
642: .   func    - a pointer to the convergence test function
643: .   ctx     - [optional] context for private data for the convergence routine
644:               (may be NULL)
645: -   destroy - [optional] destructor for the context (may be NULL;
646:               PETSC_NULL_FUNCTION in Fortran)

648:     Calling Sequence of func:
649: $   func(NEP nep,PetscInt it,PetscReal xnorm,PetscReal snorm,PetscReal fnorm,NEPConvergedReason reason*,void *fctx)

651: +   nep    - the NEP context
652: .   it     - iteration number
653: .   xnorm  - norm of the current solution
654: .   snorm  - norm of the step (difference between two consecutive solutions)
655: .   fnorm  - norm of the function (residual)
656: .   reason - (output) result of the convergence test
657: -   fctx   - optional context, as set by NEPSetConvergenceTest()

659:     Level: advanced

661: .seealso: NEPSetTolerances()
662: @*/
663: PetscErrorCode NEPSetConvergenceTest(NEP nep,PetscErrorCode (*func)(NEP,PetscInt,PetscReal,PetscReal,PetscReal,NEPConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
664: {

669:   if (nep->convergeddestroy) {
670:     (*nep->convergeddestroy)(nep->convergedctx);
671:   }
672:   nep->converged        = func;
673:   nep->convergeddestroy = destroy;
674:   nep->convergedctx     = ctx;
675:   return(0);
676: }

680: /*@
681:    NEPSetTrackAll - Specifies if the solver must compute the residual of all
682:    approximate eigenpairs or not.

684:    Logically Collective on NEP

686:    Input Parameters:
687: +  nep      - the eigensolver context
688: -  trackall - whether compute all residuals or not

690:    Notes:
691:    If the user sets trackall=PETSC_TRUE then the solver explicitly computes
692:    the residual for each eigenpair approximation. Computing the residual is
693:    usually an expensive operation and solvers commonly compute the associated
694:    residual to the first unconverged eigenpair.
695:    The options '-nep_monitor_all' and '-nep_monitor_lg_all' automatically
696:    activate this option.

698:    Level: intermediate

700: .seealso: NEPGetTrackAll()
701: @*/
702: PetscErrorCode NEPSetTrackAll(NEP nep,PetscBool trackall)
703: {
707:   nep->trackall = trackall;
708:   return(0);
709: }

713: /*@
714:    NEPGetTrackAll - Returns the flag indicating whether all residual norms must
715:    be computed or not.

717:    Not Collective

719:    Input Parameter:
720: .  nep - the eigensolver context

722:    Output Parameter:
723: .  trackall - the returned flag

725:    Level: intermediate

727: .seealso: NEPSetTrackAll()
728: @*/
729: PetscErrorCode NEPGetTrackAll(NEP nep,PetscBool *trackall)
730: {
734:   *trackall = nep->trackall;
735:   return(0);
736: }

740: /*@
741:    NEPSetRefine - Specifies the refinement type (and options) to be used
742:    after the solve.

744:    Logically Collective on NEP

746:    Input Parameters:
747: +  nep    - the nonlinear eigensolver context
748: .  refine - refinement type
749: .  tol    - the convergence tolerance
750: -  its    - maximum number of refinement iterations

752:    Options Database Keys:
753: +  -nep_refine <type> - refinement type, one of <none,simple,multiple>
754: .  -nep_refine_tol <tol> - the tolerance
755: -  -nep_refine_its <its> - number of iterations

757:    Notes:
758:    By default, iterative refinement is disabled, since it may be very
759:    costly. There are two possible refinement strategies: simple and multiple.
760:    The simple approach performs iterative refinement on each of the
761:    converged eigenpairs individually, whereas the multiple strategy works
762:    with the invariant pair as a whole, refining all eigenpairs simultaneously.
763:    The latter may be required for the case of multiple eigenvalues.

765:    The tol and its parameters specify the stopping criterion. In the simple
766:    method, refinement continues until the residual of each eigenpair is
767:    below the tolerance (tol defaults to the NEP tol, but may be set to a
768:    different value). In contrast, the multiple method simply performs its
769:    refinement iterations (just one by default).

771:    Level: intermediate

773: .seealso: NEPGetRefine()
774: @*/
775: PetscErrorCode NEPSetRefine(NEP nep,NEPRefine refine,PetscReal tol,PetscInt its)
776: {
782:   nep->refine = refine;
783:   if (refine) {  /* process parameters only if not REFINE_NONE */
784:     if (tol == PETSC_DEFAULT || tol == PETSC_DECIDE) {
785:       nep->reftol = nep->rtol;
786:     } else {
787:       if (tol<=0.0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
788:       nep->reftol = tol;
789:     }
790:     if (its==PETSC_DECIDE || its==PETSC_DEFAULT) {
791:       nep->rits = PETSC_DEFAULT;
792:     } else {
793:       if (its<=0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be > 0");
794:       nep->rits = its;
795:     }
796:   }
797:   nep->state = NEP_STATE_INITIAL;
798:   return(0);
799: }

803: /*@
804:    NEPGetRefine - Gets the refinement strategy used by the NEP object, and the
805:    associated parameters.

807:    Not Collective

809:    Input Parameter:
810: .  nep - the nonlinear eigensolver context

812:    Output Parameters:
813: +  refine - refinement type
814: .  tol    - the convergence tolerance
815: -  its    - maximum number of refinement iterations

817:    Level: intermediate

819:    Note:
820:    The user can specify NULL for any parameter that is not needed.

822: .seealso: NEPSetRefine()
823: @*/
824: PetscErrorCode NEPGetRefine(NEP nep,NEPRefine *refine,PetscReal *tol,PetscInt *its)
825: {
828:   if (refine) *refine = nep->refine;
829:   if (tol)    *tol    = nep->reftol;
830:   if (its)    *its    = nep->rits;
831:   return(0);
832: }

836: /*@C
837:    NEPSetOptionsPrefix - Sets the prefix used for searching for all
838:    NEP options in the database.

840:    Logically Collective on NEP

842:    Input Parameters:
843: +  nep - the nonlinear eigensolver context
844: -  prefix - the prefix string to prepend to all NEP option requests

846:    Notes:
847:    A hyphen (-) must NOT be given at the beginning of the prefix name.
848:    The first character of all runtime options is AUTOMATICALLY the
849:    hyphen.

851:    For example, to distinguish between the runtime options for two
852:    different NEP contexts, one could call
853: .vb
854:       NEPSetOptionsPrefix(nep1,"neig1_")
855:       NEPSetOptionsPrefix(nep2,"neig2_")
856: .ve

858:    Level: advanced

860: .seealso: NEPAppendOptionsPrefix(), NEPGetOptionsPrefix()
861: @*/
862: PetscErrorCode NEPSetOptionsPrefix(NEP nep,const char *prefix)
863: {

868:   if (!nep->V) { NEPGetBV(nep,&nep->V); }
869:   BVSetOptionsPrefix(nep->V,prefix);
870:   if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
871:   DSSetOptionsPrefix(nep->ds,prefix);
872:   if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
873:   RGSetOptionsPrefix(nep->rg,prefix);
874:   if (!nep->ksp) { NEPGetKSP(nep,&nep->ksp); }
875:   KSPSetOptionsPrefix(nep->ksp,prefix);
876:   KSPAppendOptionsPrefix(nep->ksp,"nep_");
877:   PetscObjectSetOptionsPrefix((PetscObject)nep,prefix);
878:   return(0);
879: }

883: /*@C
884:    NEPAppendOptionsPrefix - Appends to the prefix used for searching for all
885:    NEP options in the database.

887:    Logically Collective on NEP

889:    Input Parameters:
890: +  nep - the nonlinear eigensolver context
891: -  prefix - the prefix string to prepend to all NEP option requests

893:    Notes:
894:    A hyphen (-) must NOT be given at the beginning of the prefix name.
895:    The first character of all runtime options is AUTOMATICALLY the hyphen.

897:    Level: advanced

899: .seealso: NEPSetOptionsPrefix(), NEPGetOptionsPrefix()
900: @*/
901: PetscErrorCode NEPAppendOptionsPrefix(NEP nep,const char *prefix)
902: {

907:   if (!nep->V) { NEPGetBV(nep,&nep->V); }
908:   BVSetOptionsPrefix(nep->V,prefix);
909:   if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
910:   DSSetOptionsPrefix(nep->ds,prefix);
911:   if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
912:   RGSetOptionsPrefix(nep->rg,prefix);
913:   if (!nep->ksp) { NEPGetKSP(nep,&nep->ksp); }
914:   KSPSetOptionsPrefix(nep->ksp,prefix);
915:   KSPAppendOptionsPrefix(nep->ksp,"nep_");
916:   PetscObjectAppendOptionsPrefix((PetscObject)nep,prefix);
917:   return(0);
918: }

922: /*@C
923:    NEPGetOptionsPrefix - Gets the prefix used for searching for all
924:    NEP options in the database.

926:    Not Collective

928:    Input Parameters:
929: .  nep - the nonlinear eigensolver context

931:    Output Parameters:
932: .  prefix - pointer to the prefix string used is returned

934:    Notes: On the fortran side, the user should pass in a string 'prefix' of
935:    sufficient length to hold the prefix.

937:    Level: advanced

939: .seealso: NEPSetOptionsPrefix(), NEPAppendOptionsPrefix()
940: @*/
941: PetscErrorCode NEPGetOptionsPrefix(NEP nep,const char *prefix[])
942: {

948:   PetscObjectGetOptionsPrefix((PetscObject)nep,prefix);
949:   return(0);
950: }