Actual source code: epsopts.c

slepc-3.16.1 2021-11-17
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2021, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    EPS routines related to options that can be set via the command-line
 12:    or procedurally.
 13: */

 15: #include <slepc/private/epsimpl.h>
 16: #include <petscdraw.h>

 18: /*@C
 19:    EPSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
 20:    indicated by the user.

 22:    Collective on eps

 24:    Input Parameters:
 25: +  eps      - the eigensolver context
 26: .  opt      - the command line option for this monitor
 27: .  name     - the monitor type one is seeking
 28: .  ctx      - an optional user context for the monitor, or NULL
 29: -  trackall - whether this monitor tracks all eigenvalues or not

 31:    Level: developer

 33: .seealso: EPSMonitorSet(), EPSSetTrackAll()
 34: @*/
 35: PetscErrorCode EPSMonitorSetFromOptions(EPS eps,const char opt[],const char name[],void *ctx,PetscBool trackall)
 36: {
 37:   PetscErrorCode       (*mfunc)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
 38:   PetscErrorCode       (*cfunc)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**);
 39:   PetscErrorCode       (*dfunc)(PetscViewerAndFormat**);
 40:   PetscViewerAndFormat *vf;
 41:   PetscViewer          viewer;
 42:   PetscViewerFormat    format;
 43:   PetscViewerType      vtype;
 44:   char                 key[PETSC_MAX_PATH_LEN];
 45:   PetscBool            flg;
 46:   PetscErrorCode       ierr;

 49:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,opt,&viewer,&format,&flg);
 50:   if (!flg) return(0);

 52:   PetscViewerGetType(viewer,&vtype);
 53:   SlepcMonitorMakeKey_Internal(name,vtype,format,key);
 54:   PetscFunctionListFind(EPSMonitorList,key,&mfunc);
 55:   if (!mfunc) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Specified viewer and format not supported");
 56:   PetscFunctionListFind(EPSMonitorCreateList,key,&cfunc);
 57:   PetscFunctionListFind(EPSMonitorDestroyList,key,&dfunc);
 58:   if (!cfunc) cfunc = PetscViewerAndFormatCreate_Internal;
 59:   if (!dfunc) dfunc = PetscViewerAndFormatDestroy;

 61:   (*cfunc)(viewer,format,ctx,&vf);
 62:   PetscObjectDereference((PetscObject)viewer);
 63:   EPSMonitorSet(eps,mfunc,vf,(PetscErrorCode(*)(void **))dfunc);
 64:   if (trackall) {
 65:     EPSSetTrackAll(eps,PETSC_TRUE);
 66:   }
 67:   return(0);
 68: }

 70: /*@
 71:    EPSSetFromOptions - Sets EPS options from the options database.
 72:    This routine must be called before EPSSetUp() if the user is to be
 73:    allowed to set the solver type.

 75:    Collective on eps

 77:    Input Parameters:
 78: .  eps - the eigensolver context

 80:    Notes:
 81:    To see all options, run your program with the -help option.

 83:    Level: beginner
 84: @*/
 85: PetscErrorCode EPSSetFromOptions(EPS eps)
 86: {
 88:   char           type[256];
 89:   PetscBool      set,flg,flg1,flg2,flg3,bval;
 90:   PetscReal      r,array[2]={0,0};
 91:   PetscScalar    s;
 92:   PetscInt       i,j,k;
 93:   EPSBalance     bal;

 97:   EPSRegisterAll();
 98:   PetscObjectOptionsBegin((PetscObject)eps);
 99:     PetscOptionsFList("-eps_type","Eigensolver method","EPSSetType",EPSList,(char*)(((PetscObject)eps)->type_name?((PetscObject)eps)->type_name:EPSKRYLOVSCHUR),type,sizeof(type),&flg);
100:     if (flg) {
101:       EPSSetType(eps,type);
102:     } else if (!((PetscObject)eps)->type_name) {
103:       EPSSetType(eps,EPSKRYLOVSCHUR);
104:     }

106:     PetscOptionsBoolGroupBegin("-eps_hermitian","Hermitian eigenvalue problem","EPSSetProblemType",&flg);
107:     if (flg) { EPSSetProblemType(eps,EPS_HEP); }
108:     PetscOptionsBoolGroup("-eps_gen_hermitian","Generalized Hermitian eigenvalue problem","EPSSetProblemType",&flg);
109:     if (flg) { EPSSetProblemType(eps,EPS_GHEP); }
110:     PetscOptionsBoolGroup("-eps_non_hermitian","Non-Hermitian eigenvalue problem","EPSSetProblemType",&flg);
111:     if (flg) { EPSSetProblemType(eps,EPS_NHEP); }
112:     PetscOptionsBoolGroup("-eps_gen_non_hermitian","Generalized non-Hermitian eigenvalue problem","EPSSetProblemType",&flg);
113:     if (flg) { EPSSetProblemType(eps,EPS_GNHEP); }
114:     PetscOptionsBoolGroup("-eps_pos_gen_non_hermitian","Generalized non-Hermitian eigenvalue problem with positive semi-definite B","EPSSetProblemType",&flg);
115:     if (flg) { EPSSetProblemType(eps,EPS_PGNHEP); }
116:     PetscOptionsBoolGroupEnd("-eps_gen_indefinite","Generalized Hermitian-indefinite eigenvalue problem","EPSSetProblemType",&flg);
117:     if (flg) { EPSSetProblemType(eps,EPS_GHIEP); }

119:     PetscOptionsBoolGroupBegin("-eps_ritz","Rayleigh-Ritz extraction","EPSSetExtraction",&flg);
120:     if (flg) { EPSSetExtraction(eps,EPS_RITZ); }
121:     PetscOptionsBoolGroup("-eps_harmonic","Harmonic Ritz extraction","EPSSetExtraction",&flg);
122:     if (flg) { EPSSetExtraction(eps,EPS_HARMONIC); }
123:     PetscOptionsBoolGroup("-eps_harmonic_relative","Relative harmonic Ritz extraction","EPSSetExtraction",&flg);
124:     if (flg) { EPSSetExtraction(eps,EPS_HARMONIC_RELATIVE); }
125:     PetscOptionsBoolGroup("-eps_harmonic_right","Right harmonic Ritz extraction","EPSSetExtraction",&flg);
126:     if (flg) { EPSSetExtraction(eps,EPS_HARMONIC_RIGHT); }
127:     PetscOptionsBoolGroup("-eps_harmonic_largest","Largest harmonic Ritz extraction","EPSSetExtraction",&flg);
128:     if (flg) { EPSSetExtraction(eps,EPS_HARMONIC_LARGEST); }
129:     PetscOptionsBoolGroup("-eps_refined","Refined Ritz extraction","EPSSetExtraction",&flg);
130:     if (flg) { EPSSetExtraction(eps,EPS_REFINED); }
131:     PetscOptionsBoolGroupEnd("-eps_refined_harmonic","Refined harmonic Ritz extraction","EPSSetExtraction",&flg);
132:     if (flg) { EPSSetExtraction(eps,EPS_REFINED_HARMONIC); }

134:     bal = eps->balance;
135:     PetscOptionsEnum("-eps_balance","Balancing method","EPSSetBalance",EPSBalanceTypes,(PetscEnum)bal,(PetscEnum*)&bal,&flg1);
136:     j = eps->balance_its;
137:     PetscOptionsInt("-eps_balance_its","Number of iterations in balancing","EPSSetBalance",eps->balance_its,&j,&flg2);
138:     r = eps->balance_cutoff;
139:     PetscOptionsReal("-eps_balance_cutoff","Cutoff value in balancing","EPSSetBalance",eps->balance_cutoff,&r,&flg3);
140:     if (flg1 || flg2 || flg3) { EPSSetBalance(eps,bal,j,r); }

142:     i = eps->max_it;
143:     PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&i,&flg1);
144:     r = eps->tol;
145:     PetscOptionsReal("-eps_tol","Tolerance","EPSSetTolerances",SlepcDefaultTol(eps->tol),&r,&flg2);
146:     if (flg1 || flg2) { EPSSetTolerances(eps,r,i); }

148:     PetscOptionsBoolGroupBegin("-eps_conv_rel","Relative error convergence test","EPSSetConvergenceTest",&flg);
149:     if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_REL); }
150:     PetscOptionsBoolGroup("-eps_conv_norm","Convergence test relative to the eigenvalue and the matrix norms","EPSSetConvergenceTest",&flg);
151:     if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_NORM); }
152:     PetscOptionsBoolGroup("-eps_conv_abs","Absolute error convergence test","EPSSetConvergenceTest",&flg);
153:     if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_ABS); }
154:     PetscOptionsBoolGroupEnd("-eps_conv_user","User-defined convergence test","EPSSetConvergenceTest",&flg);
155:     if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_USER); }

157:     PetscOptionsBoolGroupBegin("-eps_stop_basic","Stop iteration if all eigenvalues converged or max_it reached","EPSSetStoppingTest",&flg);
158:     if (flg) { EPSSetStoppingTest(eps,EPS_STOP_BASIC); }
159:     PetscOptionsBoolGroupEnd("-eps_stop_user","User-defined stopping test","EPSSetStoppingTest",&flg);
160:     if (flg) { EPSSetStoppingTest(eps,EPS_STOP_USER); }

162:     i = eps->nev;
163:     PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&i,&flg1);
164:     j = eps->ncv;
165:     PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&j,&flg2);
166:     k = eps->mpd;
167:     PetscOptionsInt("-eps_mpd","Maximum dimension of projected problem","EPSSetDimensions",eps->mpd,&k,&flg3);
168:     if (flg1 || flg2 || flg3) { EPSSetDimensions(eps,i,j,k); }

170:     PetscOptionsBoolGroupBegin("-eps_largest_magnitude","Compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
171:     if (flg) { EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE); }
172:     PetscOptionsBoolGroup("-eps_smallest_magnitude","Compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
173:     if (flg) { EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE); }
174:     PetscOptionsBoolGroup("-eps_largest_real","Compute eigenvalues with largest real parts","EPSSetWhichEigenpairs",&flg);
175:     if (flg) { EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL); }
176:     PetscOptionsBoolGroup("-eps_smallest_real","Compute eigenvalues with smallest real parts","EPSSetWhichEigenpairs",&flg);
177:     if (flg) { EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL); }
178:     PetscOptionsBoolGroup("-eps_largest_imaginary","Compute eigenvalues with largest imaginary parts","EPSSetWhichEigenpairs",&flg);
179:     if (flg) { EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY); }
180:     PetscOptionsBoolGroup("-eps_smallest_imaginary","Compute eigenvalues with smallest imaginary parts","EPSSetWhichEigenpairs",&flg);
181:     if (flg) { EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY); }
182:     PetscOptionsBoolGroup("-eps_target_magnitude","Compute eigenvalues closest to target","EPSSetWhichEigenpairs",&flg);
183:     if (flg) { EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE); }
184:     PetscOptionsBoolGroup("-eps_target_real","Compute eigenvalues with real parts closest to target","EPSSetWhichEigenpairs",&flg);
185:     if (flg) { EPSSetWhichEigenpairs(eps,EPS_TARGET_REAL); }
186:     PetscOptionsBoolGroup("-eps_target_imaginary","Compute eigenvalues with imaginary parts closest to target","EPSSetWhichEigenpairs",&flg);
187:     if (flg) { EPSSetWhichEigenpairs(eps,EPS_TARGET_IMAGINARY); }
188:     PetscOptionsBoolGroupEnd("-eps_all","Compute all eigenvalues in an interval or a region","EPSSetWhichEigenpairs",&flg);
189:     if (flg) { EPSSetWhichEigenpairs(eps,EPS_ALL); }

191:     PetscOptionsScalar("-eps_target","Value of the target","EPSSetTarget",eps->target,&s,&flg);
192:     if (flg) {
193:       if (eps->which!=EPS_TARGET_REAL && eps->which!=EPS_TARGET_IMAGINARY) {
194:         EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);
195:       }
196:       EPSSetTarget(eps,s);
197:     }

199:     k = 2;
200:     PetscOptionsRealArray("-eps_interval","Computational interval (two real values separated with a comma without spaces)","EPSSetInterval",array,&k,&flg);
201:     if (flg) {
202:       if (k<2) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_SIZ,"Must pass two values in -eps_interval (comma-separated without spaces)");
203:       EPSSetWhichEigenpairs(eps,EPS_ALL);
204:       EPSSetInterval(eps,array[0],array[1]);
205:     }

207:     PetscOptionsBool("-eps_true_residual","Compute true residuals explicitly","EPSSetTrueResidual",eps->trueres,&eps->trueres,NULL);
208:     PetscOptionsBool("-eps_purify","Postprocess eigenvectors for purification","EPSSetPurify",eps->purify,&bval,&flg);
209:     if (flg) { EPSSetPurify(eps,bval); }
210:     PetscOptionsBool("-eps_two_sided","Use two-sided variant (to compute left eigenvectors)","EPSSetTwoSided",eps->twosided,&bval,&flg);
211:     if (flg) { EPSSetTwoSided(eps,bval); }

213:     /* -----------------------------------------------------------------------*/
214:     /*
215:       Cancels all monitors hardwired into code before call to EPSSetFromOptions()
216:     */
217:     PetscOptionsBool("-eps_monitor_cancel","Remove any hardwired monitor routines","EPSMonitorCancel",PETSC_FALSE,&flg,&set);
218:     if (set && flg) { EPSMonitorCancel(eps); }
219:     EPSMonitorSetFromOptions(eps,"-eps_monitor","first_approximation",NULL,PETSC_FALSE);
220:     EPSMonitorSetFromOptions(eps,"-eps_monitor_all","all_approximations",NULL,PETSC_TRUE);
221:     EPSMonitorSetFromOptions(eps,"-eps_monitor_conv","convergence_history",NULL,PETSC_FALSE);

223:     /* -----------------------------------------------------------------------*/
224:     PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",NULL);
225:     PetscOptionsName("-eps_view_vectors","View computed eigenvectors","EPSVectorsView",NULL);
226:     PetscOptionsName("-eps_view_values","View computed eigenvalues","EPSValuesView",NULL);
227:     PetscOptionsName("-eps_converged_reason","Print reason for convergence, and number of iterations","EPSConvergedReasonView",NULL);
228:     PetscOptionsName("-eps_error_absolute","Print absolute errors of each eigenpair","EPSErrorView",NULL);
229:     PetscOptionsName("-eps_error_relative","Print relative errors of each eigenpair","EPSErrorView",NULL);
230:     PetscOptionsName("-eps_error_backward","Print backward errors of each eigenpair","EPSErrorView",NULL);

232:     if (eps->ops->setfromoptions) {
233:       (*eps->ops->setfromoptions)(PetscOptionsObject,eps);
234:     }
235:     PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)eps);
236:   PetscOptionsEnd();

238:   if (!eps->V) { EPSGetBV(eps,&eps->V); }
239:   BVSetFromOptions(eps->V);
240:   if (!eps->rg) { EPSGetRG(eps,&eps->rg); }
241:   RGSetFromOptions(eps->rg);
242:   if (eps->useds) {
243:     if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
244:     DSSetFromOptions(eps->ds);
245:   }
246:   if (!eps->st) { EPSGetST(eps,&eps->st); }
247:   EPSSetDefaultST(eps);
248:   STSetFromOptions(eps->st);
249:   return(0);
250: }

252: /*@C
253:    EPSGetTolerances - Gets the tolerance and maximum iteration count used
254:    by the EPS convergence tests.

256:    Not Collective

258:    Input Parameter:
259: .  eps - the eigensolver context

261:    Output Parameters:
262: +  tol - the convergence tolerance
263: -  maxits - maximum number of iterations

265:    Notes:
266:    The user can specify NULL for any parameter that is not needed.

268:    Level: intermediate

270: .seealso: EPSSetTolerances()
271: @*/
272: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,PetscInt *maxits)
273: {
276:   if (tol)    *tol    = eps->tol;
277:   if (maxits) *maxits = eps->max_it;
278:   return(0);
279: }

281: /*@
282:    EPSSetTolerances - Sets the tolerance and maximum iteration count used
283:    by the EPS convergence tests.

285:    Logically Collective on eps

287:    Input Parameters:
288: +  eps - the eigensolver context
289: .  tol - the convergence tolerance
290: -  maxits - maximum number of iterations to use

292:    Options Database Keys:
293: +  -eps_tol <tol> - Sets the convergence tolerance
294: -  -eps_max_it <maxits> - Sets the maximum number of iterations allowed

296:    Notes:
297:    Use PETSC_DEFAULT for either argument to assign a reasonably good value.

299:    Level: intermediate

301: .seealso: EPSGetTolerances()
302: @*/
303: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,PetscInt maxits)
304: {
309:   if (tol == PETSC_DEFAULT) {
310:     eps->tol   = PETSC_DEFAULT;
311:     eps->state = EPS_STATE_INITIAL;
312:   } else {
313:     if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
314:     eps->tol = tol;
315:   }
316:   if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
317:     eps->max_it = PETSC_DEFAULT;
318:     eps->state  = EPS_STATE_INITIAL;
319:   } else {
320:     if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
321:     eps->max_it = maxits;
322:   }
323:   return(0);
324: }

326: /*@C
327:    EPSGetDimensions - Gets the number of eigenvalues to compute
328:    and the dimension of the subspace.

330:    Not Collective

332:    Input Parameter:
333: .  eps - the eigensolver context

335:    Output Parameters:
336: +  nev - number of eigenvalues to compute
337: .  ncv - the maximum dimension of the subspace to be used by the solver
338: -  mpd - the maximum dimension allowed for the projected problem

340:    Level: intermediate

342: .seealso: EPSSetDimensions()
343: @*/
344: PetscErrorCode EPSGetDimensions(EPS eps,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
345: {
348:   if (nev) *nev = eps->nev;
349:   if (ncv) *ncv = eps->ncv;
350:   if (mpd) *mpd = eps->mpd;
351:   return(0);
352: }

354: /*@
355:    EPSSetDimensions - Sets the number of eigenvalues to compute
356:    and the dimension of the subspace.

358:    Logically Collective on eps

360:    Input Parameters:
361: +  eps - the eigensolver context
362: .  nev - number of eigenvalues to compute
363: .  ncv - the maximum dimension of the subspace to be used by the solver
364: -  mpd - the maximum dimension allowed for the projected problem

366:    Options Database Keys:
367: +  -eps_nev <nev> - Sets the number of eigenvalues
368: .  -eps_ncv <ncv> - Sets the dimension of the subspace
369: -  -eps_mpd <mpd> - Sets the maximum projected dimension

371:    Notes:
372:    Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
373:    dependent on the solution method.

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

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

384:    When computing all eigenvalues in an interval, see EPSSetInterval(), these
385:    parameters lose relevance, and tuning must be done with
386:    EPSKrylovSchurSetDimensions().

388:    Level: intermediate

390: .seealso: EPSGetDimensions(), EPSSetInterval(), EPSKrylovSchurSetDimensions()
391: @*/
392: PetscErrorCode EPSSetDimensions(EPS eps,PetscInt nev,PetscInt ncv,PetscInt mpd)
393: {
399:   if (nev<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
400:   eps->nev = nev;
401:   if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
402:     eps->ncv = PETSC_DEFAULT;
403:   } else {
404:     if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
405:     eps->ncv = ncv;
406:   }
407:   if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
408:     eps->mpd = PETSC_DEFAULT;
409:   } else {
410:     if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
411:     eps->mpd = mpd;
412:   }
413:   eps->state = EPS_STATE_INITIAL;
414:   return(0);
415: }

417: /*@
418:    EPSSetWhichEigenpairs - Specifies which portion of the spectrum is
419:    to be sought.

421:    Logically Collective on eps

423:    Input Parameters:
424: +  eps   - eigensolver context obtained from EPSCreate()
425: -  which - the portion of the spectrum to be sought

427:    Possible values:
428:    The parameter 'which' can have one of these values

430: +     EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
431: .     EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
432: .     EPS_LARGEST_REAL - largest real parts
433: .     EPS_SMALLEST_REAL - smallest real parts
434: .     EPS_LARGEST_IMAGINARY - largest imaginary parts
435: .     EPS_SMALLEST_IMAGINARY - smallest imaginary parts
436: .     EPS_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
437: .     EPS_TARGET_REAL - eigenvalues with real part closest to target
438: .     EPS_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
439: .     EPS_ALL - all eigenvalues contained in a given interval or region
440: -     EPS_WHICH_USER - user defined ordering set with EPSSetEigenvalueComparison()

442:    Options Database Keys:
443: +   -eps_largest_magnitude - Sets largest eigenvalues in magnitude
444: .   -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
445: .   -eps_largest_real - Sets largest real parts
446: .   -eps_smallest_real - Sets smallest real parts
447: .   -eps_largest_imaginary - Sets largest imaginary parts
448: .   -eps_smallest_imaginary - Sets smallest imaginary parts
449: .   -eps_target_magnitude - Sets eigenvalues closest to target
450: .   -eps_target_real - Sets real parts closest to target
451: .   -eps_target_imaginary - Sets imaginary parts closest to target
452: -   -eps_all - Sets all eigenvalues in an interval or region

454:    Notes:
455:    Not all eigensolvers implemented in EPS account for all the possible values
456:    stated above. Also, some values make sense only for certain types of
457:    problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY
458:    and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part
459:    for eigenvalue selection.

461:    The target is a scalar value provided with EPSSetTarget().

463:    The criterion EPS_TARGET_IMAGINARY is available only in case PETSc and
464:    SLEPc have been built with complex scalars.

466:    EPS_ALL is intended for use in combination with an interval (see
467:    EPSSetInterval()), when all eigenvalues within the interval are requested,
468:    or in the context of the CISS solver for computing all eigenvalues in a region.
469:    In those cases, the number of eigenvalues is unknown, so the nev parameter
470:    has a different sense, see EPSSetDimensions().

472:    Level: intermediate

474: .seealso: EPSGetWhichEigenpairs(), EPSSetTarget(), EPSSetInterval(),
475:           EPSSetDimensions(), EPSSetEigenvalueComparison(), EPSWhich
476: @*/
477: PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)
478: {
482:   switch (which) {
483:     case EPS_LARGEST_MAGNITUDE:
484:     case EPS_SMALLEST_MAGNITUDE:
485:     case EPS_LARGEST_REAL:
486:     case EPS_SMALLEST_REAL:
487:     case EPS_LARGEST_IMAGINARY:
488:     case EPS_SMALLEST_IMAGINARY:
489:     case EPS_TARGET_MAGNITUDE:
490:     case EPS_TARGET_REAL:
491: #if defined(PETSC_USE_COMPLEX)
492:     case EPS_TARGET_IMAGINARY:
493: #endif
494:     case EPS_ALL:
495:     case EPS_WHICH_USER:
496:       if (eps->which != which) {
497:         eps->state = EPS_STATE_INITIAL;
498:         eps->which = which;
499:       }
500:       break;
501: #if !defined(PETSC_USE_COMPLEX)
502:     case EPS_TARGET_IMAGINARY:
503:       SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"EPS_TARGET_IMAGINARY can be used only with complex scalars");
504: #endif
505:     default:
506:       SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
507:   }
508:   return(0);
509: }

511: /*@
512:    EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be
513:    sought.

515:    Not Collective

517:    Input Parameter:
518: .  eps - eigensolver context obtained from EPSCreate()

520:    Output Parameter:
521: .  which - the portion of the spectrum to be sought

523:    Notes:
524:    See EPSSetWhichEigenpairs() for possible values of 'which'.

526:    Level: intermediate

528: .seealso: EPSSetWhichEigenpairs(), EPSWhich
529: @*/
530: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)
531: {
535:   *which = eps->which;
536:   return(0);
537: }

539: /*@C
540:    EPSSetEigenvalueComparison - Specifies the eigenvalue comparison function
541:    when EPSSetWhichEigenpairs() is set to EPS_WHICH_USER.

543:    Logically Collective on eps

545:    Input Parameters:
546: +  eps  - eigensolver context obtained from EPSCreate()
547: .  func - a pointer to the comparison function
548: -  ctx  - a context pointer (the last parameter to the comparison function)

550:    Calling Sequence of func:
551: $   func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)

553: +   ar     - real part of the 1st eigenvalue
554: .   ai     - imaginary part of the 1st eigenvalue
555: .   br     - real part of the 2nd eigenvalue
556: .   bi     - imaginary part of the 2nd eigenvalue
557: .   res    - result of comparison
558: -   ctx    - optional context, as set by EPSSetEigenvalueComparison()

560:    Note:
561:    The returning parameter 'res' can be
562: +  negative - if the 1st eigenvalue is preferred to the 2st one
563: .  zero     - if both eigenvalues are equally preferred
564: -  positive - if the 2st eigenvalue is preferred to the 1st one

566:    Level: advanced

568: .seealso: EPSSetWhichEigenpairs(), EPSWhich
569: @*/
570: PetscErrorCode EPSSetEigenvalueComparison(EPS eps,PetscErrorCode (*func)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)
571: {
574:   eps->sc->comparison    = func;
575:   eps->sc->comparisonctx = ctx;
576:   eps->which             = EPS_WHICH_USER;
577:   return(0);
578: }

580: /*@C
581:    EPSSetArbitrarySelection - Specifies a function intended to look for
582:    eigenvalues according to an arbitrary selection criterion. This criterion
583:    can be based on a computation involving the current eigenvector approximation.

585:    Logically Collective on eps

587:    Input Parameters:
588: +  eps  - eigensolver context obtained from EPSCreate()
589: .  func - a pointer to the evaluation function
590: -  ctx  - a context pointer (the last parameter to the evaluation function)

592:    Calling Sequence of func:
593: $   func(PetscScalar er,PetscScalar ei,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx)

595: +   er     - real part of the current eigenvalue approximation
596: .   ei     - imaginary part of the current eigenvalue approximation
597: .   xr     - real part of the current eigenvector approximation
598: .   xi     - imaginary part of the current eigenvector approximation
599: .   rr     - result of evaluation (real part)
600: .   ri     - result of evaluation (imaginary part)
601: -   ctx    - optional context, as set by EPSSetArbitrarySelection()

603:    Notes:
604:    This provides a mechanism to select eigenpairs by evaluating a user-defined
605:    function. When a function has been provided, the default selection based on
606:    sorting the eigenvalues is replaced by the sorting of the results of this
607:    function (with the same sorting criterion given in EPSSetWhichEigenpairs()).

609:    For instance, suppose you want to compute those eigenvectors that maximize
610:    a certain computable expression. Then implement the computation using
611:    the arguments xr and xi, and return the result in rr. Then set the standard
612:    sorting by magnitude so that the eigenpair with largest value of rr is
613:    selected.

615:    This evaluation function is collective, that is, all processes call it and
616:    it can use collective operations; furthermore, the computed result must
617:    be the same in all processes.

619:    The result of func is expressed as a complex number so that it is possible to
620:    use the standard eigenvalue sorting functions, but normally only rr is used.
621:    Set ri to zero unless it is meaningful in your application.

623:    Level: advanced

625: .seealso: EPSSetWhichEigenpairs()
626: @*/
627: PetscErrorCode EPSSetArbitrarySelection(EPS eps,PetscErrorCode (*func)(PetscScalar,PetscScalar,Vec,Vec,PetscScalar*,PetscScalar*,void*),void* ctx)
628: {
631:   eps->arbitrary    = func;
632:   eps->arbitraryctx = ctx;
633:   eps->state        = EPS_STATE_INITIAL;
634:   return(0);
635: }

637: /*@C
638:    EPSSetConvergenceTestFunction - Sets a function to compute the error estimate
639:    used in the convergence test.

641:    Logically Collective on eps

643:    Input Parameters:
644: +  eps     - eigensolver context obtained from EPSCreate()
645: .  func    - a pointer to the convergence test function
646: .  ctx     - context for private data for the convergence routine (may be null)
647: -  destroy - a routine for destroying the context (may be null)

649:    Calling Sequence of func:
650: $   func(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)

652: +   eps    - eigensolver context obtained from EPSCreate()
653: .   eigr   - real part of the eigenvalue
654: .   eigi   - imaginary part of the eigenvalue
655: .   res    - residual norm associated to the eigenpair
656: .   errest - (output) computed error estimate
657: -   ctx    - optional context, as set by EPSSetConvergenceTestFunction()

659:    Note:
660:    If the error estimate returned by the convergence test function is less than
661:    the tolerance, then the eigenvalue is accepted as converged.

663:    Level: advanced

665: .seealso: EPSSetConvergenceTest(), EPSSetTolerances()
666: @*/
667: PetscErrorCode EPSSetConvergenceTestFunction(EPS eps,PetscErrorCode (*func)(EPS,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
668: {

673:   if (eps->convergeddestroy) {
674:     (*eps->convergeddestroy)(eps->convergedctx);
675:   }
676:   eps->convergeduser    = func;
677:   eps->convergeddestroy = destroy;
678:   eps->convergedctx     = ctx;
679:   if (func == EPSConvergedRelative) eps->conv = EPS_CONV_REL;
680:   else if (func == EPSConvergedNorm) eps->conv = EPS_CONV_NORM;
681:   else if (func == EPSConvergedAbsolute) eps->conv = EPS_CONV_ABS;
682:   else {
683:     eps->conv      = EPS_CONV_USER;
684:     eps->converged = eps->convergeduser;
685:   }
686:   return(0);
687: }

689: /*@
690:    EPSSetConvergenceTest - Specifies how to compute the error estimate
691:    used in the convergence test.

693:    Logically Collective on eps

695:    Input Parameters:
696: +  eps  - eigensolver context obtained from EPSCreate()
697: -  conv - the type of convergence test

699:    Options Database Keys:
700: +  -eps_conv_abs  - Sets the absolute convergence test
701: .  -eps_conv_rel  - Sets the convergence test relative to the eigenvalue
702: .  -eps_conv_norm - Sets the convergence test relative to the matrix norms
703: -  -eps_conv_user - Selects the user-defined convergence test

705:    Note:
706:    The parameter 'conv' can have one of these values
707: +     EPS_CONV_ABS  - absolute error ||r||
708: .     EPS_CONV_REL  - error relative to the eigenvalue l, ||r||/|l|
709: .     EPS_CONV_NORM - error relative to the matrix norms, ||r||/(||A||+|l|*||B||)
710: -     EPS_CONV_USER - function set by EPSSetConvergenceTestFunction()

712:    Level: intermediate

714: .seealso: EPSGetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSSetStoppingTest(), EPSConv
715: @*/
716: PetscErrorCode EPSSetConvergenceTest(EPS eps,EPSConv conv)
717: {
721:   switch (conv) {
722:     case EPS_CONV_ABS:  eps->converged = EPSConvergedAbsolute; break;
723:     case EPS_CONV_REL:  eps->converged = EPSConvergedRelative; break;
724:     case EPS_CONV_NORM: eps->converged = EPSConvergedNorm; break;
725:     case EPS_CONV_USER:
726:       if (!eps->convergeduser) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ORDER,"Must call EPSSetConvergenceTestFunction() first");
727:       eps->converged = eps->convergeduser;
728:       break;
729:     default:
730:       SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
731:   }
732:   eps->conv = conv;
733:   return(0);
734: }

736: /*@
737:    EPSGetConvergenceTest - Gets the method used to compute the error estimate
738:    used in the convergence test.

740:    Not Collective

742:    Input Parameters:
743: .  eps   - eigensolver context obtained from EPSCreate()

745:    Output Parameters:
746: .  conv  - the type of convergence test

748:    Level: intermediate

750: .seealso: EPSSetConvergenceTest(), EPSConv
751: @*/
752: PetscErrorCode EPSGetConvergenceTest(EPS eps,EPSConv *conv)
753: {
757:   *conv = eps->conv;
758:   return(0);
759: }

761: /*@C
762:    EPSSetStoppingTestFunction - Sets a function to decide when to stop the outer
763:    iteration of the eigensolver.

765:    Logically Collective on eps

767:    Input Parameters:
768: +  eps     - eigensolver context obtained from EPSCreate()
769: .  func    - pointer to the stopping test function
770: .  ctx     - context for private data for the stopping routine (may be null)
771: -  destroy - a routine for destroying the context (may be null)

773:    Calling Sequence of func:
774: $   func(EPS eps,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,EPSConvergedReason *reason,void *ctx)

776: +   eps    - eigensolver context obtained from EPSCreate()
777: .   its    - current number of iterations
778: .   max_it - maximum number of iterations
779: .   nconv  - number of currently converged eigenpairs
780: .   nev    - number of requested eigenpairs
781: .   reason - (output) result of the stopping test
782: -   ctx    - optional context, as set by EPSSetStoppingTestFunction()

784:    Note:
785:    Normal usage is to first call the default routine EPSStoppingBasic() and then
786:    set reason to EPS_CONVERGED_USER if some user-defined conditions have been
787:    met. To let the eigensolver continue iterating, the result must be left as
788:    EPS_CONVERGED_ITERATING.

790:    Level: advanced

792: .seealso: EPSSetStoppingTest(), EPSStoppingBasic()
793: @*/
794: PetscErrorCode EPSSetStoppingTestFunction(EPS eps,PetscErrorCode (*func)(EPS,PetscInt,PetscInt,PetscInt,PetscInt,EPSConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
795: {

800:   if (eps->stoppingdestroy) {
801:     (*eps->stoppingdestroy)(eps->stoppingctx);
802:   }
803:   eps->stoppinguser    = func;
804:   eps->stoppingdestroy = destroy;
805:   eps->stoppingctx     = ctx;
806:   if (func == EPSStoppingBasic) eps->stop = EPS_STOP_BASIC;
807:   else {
808:     eps->stop     = EPS_STOP_USER;
809:     eps->stopping = eps->stoppinguser;
810:   }
811:   return(0);
812: }

814: /*@
815:    EPSSetStoppingTest - Specifies how to decide the termination of the outer
816:    loop of the eigensolver.

818:    Logically Collective on eps

820:    Input Parameters:
821: +  eps  - eigensolver context obtained from EPSCreate()
822: -  stop - the type of stopping test

824:    Options Database Keys:
825: +  -eps_stop_basic - Sets the default stopping test
826: -  -eps_stop_user  - Selects the user-defined stopping test

828:    Note:
829:    The parameter 'stop' can have one of these values
830: +     EPS_STOP_BASIC - default stopping test
831: -     EPS_STOP_USER  - function set by EPSSetStoppingTestFunction()

833:    Level: advanced

835: .seealso: EPSGetStoppingTest(), EPSSetStoppingTestFunction(), EPSSetConvergenceTest(), EPSStop
836: @*/
837: PetscErrorCode EPSSetStoppingTest(EPS eps,EPSStop stop)
838: {
842:   switch (stop) {
843:     case EPS_STOP_BASIC: eps->stopping = EPSStoppingBasic; break;
844:     case EPS_STOP_USER:
845:       if (!eps->stoppinguser) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ORDER,"Must call EPSSetStoppingTestFunction() first");
846:       eps->stopping = eps->stoppinguser;
847:       break;
848:     default:
849:       SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
850:   }
851:   eps->stop = stop;
852:   return(0);
853: }

855: /*@
856:    EPSGetStoppingTest - Gets the method used to decide the termination of the outer
857:    loop of the eigensolver.

859:    Not Collective

861:    Input Parameters:
862: .  eps   - eigensolver context obtained from EPSCreate()

864:    Output Parameters:
865: .  stop  - the type of stopping test

867:    Level: advanced

869: .seealso: EPSSetStoppingTest(), EPSStop
870: @*/
871: PetscErrorCode EPSGetStoppingTest(EPS eps,EPSStop *stop)
872: {
876:   *stop = eps->stop;
877:   return(0);
878: }

880: /*@
881:    EPSSetProblemType - Specifies the type of the eigenvalue problem.

883:    Logically Collective on eps

885:    Input Parameters:
886: +  eps      - the eigensolver context
887: -  type     - a known type of eigenvalue problem

889:    Options Database Keys:
890: +  -eps_hermitian - Hermitian eigenvalue problem
891: .  -eps_gen_hermitian - generalized Hermitian eigenvalue problem
892: .  -eps_non_hermitian - non-Hermitian eigenvalue problem
893: .  -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
894: .  -eps_pos_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
895:    with positive semi-definite B
896: -  -eps_gen_indefinite - generalized Hermitian-indefinite eigenvalue problem

898:    Notes:
899:    Allowed values for the problem type are: Hermitian (EPS_HEP), non-Hermitian
900:    (EPS_NHEP), generalized Hermitian (EPS_GHEP), generalized non-Hermitian
901:    (EPS_GNHEP), generalized non-Hermitian with positive semi-definite B
902:    (EPS_PGNHEP), and generalized Hermitian-indefinite (EPS_GHIEP).

904:    This function must be used to instruct SLEPc to exploit symmetry. If no
905:    problem type is specified, by default a non-Hermitian problem is assumed
906:    (either standard or generalized). If the user knows that the problem is
907:    Hermitian (i.e. A=A^H) or generalized Hermitian (i.e. A=A^H, B=B^H, and
908:    B positive definite) then it is recommended to set the problem type so
909:    that eigensolver can exploit these properties.

911:    Level: intermediate

913: .seealso: EPSSetOperators(), EPSSetType(), EPSGetProblemType(), EPSProblemType
914: @*/
915: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)
916: {
920:   if (type == eps->problem_type) return(0);
921:   switch (type) {
922:     case EPS_HEP:
923:       eps->isgeneralized = PETSC_FALSE;
924:       eps->ishermitian = PETSC_TRUE;
925:       eps->ispositive = PETSC_FALSE;
926:       break;
927:     case EPS_NHEP:
928:       eps->isgeneralized = PETSC_FALSE;
929:       eps->ishermitian = PETSC_FALSE;
930:       eps->ispositive = PETSC_FALSE;
931:       break;
932:     case EPS_GHEP:
933:       eps->isgeneralized = PETSC_TRUE;
934:       eps->ishermitian = PETSC_TRUE;
935:       eps->ispositive = PETSC_TRUE;
936:       break;
937:     case EPS_GNHEP:
938:       eps->isgeneralized = PETSC_TRUE;
939:       eps->ishermitian = PETSC_FALSE;
940:       eps->ispositive = PETSC_FALSE;
941:       break;
942:     case EPS_PGNHEP:
943:       eps->isgeneralized = PETSC_TRUE;
944:       eps->ishermitian = PETSC_FALSE;
945:       eps->ispositive = PETSC_TRUE;
946:       break;
947:     case EPS_GHIEP:
948:       eps->isgeneralized = PETSC_TRUE;
949:       eps->ishermitian = PETSC_TRUE;
950:       eps->ispositive = PETSC_FALSE;
951:       break;
952:     default:
953:       SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
954:   }
955:   eps->problem_type = type;
956:   eps->state = EPS_STATE_INITIAL;
957:   return(0);
958: }

960: /*@
961:    EPSGetProblemType - Gets the problem type from the EPS object.

963:    Not Collective

965:    Input Parameter:
966: .  eps - the eigensolver context

968:    Output Parameter:
969: .  type - the problem type

971:    Level: intermediate

973: .seealso: EPSSetProblemType(), EPSProblemType
974: @*/
975: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)
976: {
980:   *type = eps->problem_type;
981:   return(0);
982: }

984: /*@
985:    EPSSetExtraction - Specifies the type of extraction technique to be employed
986:    by the eigensolver.

988:    Logically Collective on eps

990:    Input Parameters:
991: +  eps  - the eigensolver context
992: -  extr - a known type of extraction

994:    Options Database Keys:
995: +  -eps_ritz - Rayleigh-Ritz extraction
996: .  -eps_harmonic - harmonic Ritz extraction
997: .  -eps_harmonic_relative - harmonic Ritz extraction relative to the eigenvalue
998: .  -eps_harmonic_right - harmonic Ritz extraction for rightmost eigenvalues
999: .  -eps_harmonic_largest - harmonic Ritz extraction for largest magnitude
1000:    (without target)
1001: .  -eps_refined - refined Ritz extraction
1002: -  -eps_refined_harmonic - refined harmonic Ritz extraction

1004:    Notes:
1005:    Not all eigensolvers support all types of extraction. See the SLEPc
1006:    Users Manual for details.

1008:    By default, a standard Rayleigh-Ritz extraction is used. Other extractions
1009:    may be useful when computing interior eigenvalues.

1011:    Harmonic-type extractions are used in combination with a 'target'.

1013:    Level: advanced

1015: .seealso: EPSSetTarget(), EPSGetExtraction(), EPSExtraction
1016: @*/
1017: PetscErrorCode EPSSetExtraction(EPS eps,EPSExtraction extr)
1018: {
1022:   eps->extraction = extr;
1023:   return(0);
1024: }

1026: /*@
1027:    EPSGetExtraction - Gets the extraction type used by the EPS object.

1029:    Not Collective

1031:    Input Parameter:
1032: .  eps - the eigensolver context

1034:    Output Parameter:
1035: .  extr - name of extraction type

1037:    Level: advanced

1039: .seealso: EPSSetExtraction(), EPSExtraction
1040: @*/
1041: PetscErrorCode EPSGetExtraction(EPS eps,EPSExtraction *extr)
1042: {
1046:   *extr = eps->extraction;
1047:   return(0);
1048: }

1050: /*@
1051:    EPSSetBalance - Specifies the balancing technique to be employed by the
1052:    eigensolver, and some parameters associated to it.

1054:    Logically Collective on eps

1056:    Input Parameters:
1057: +  eps    - the eigensolver context
1058: .  bal    - the balancing method, one of EPS_BALANCE_NONE, EPS_BALANCE_ONESIDE,
1059:             EPS_BALANCE_TWOSIDE, or EPS_BALANCE_USER
1060: .  its    - number of iterations of the balancing algorithm
1061: -  cutoff - cutoff value

1063:    Options Database Keys:
1064: +  -eps_balance <method> - the balancing method, where <method> is one of
1065:                            'none', 'oneside', 'twoside', or 'user'
1066: .  -eps_balance_its <its> - number of iterations
1067: -  -eps_balance_cutoff <cutoff> - cutoff value

1069:    Notes:
1070:    When balancing is enabled, the solver works implicitly with matrix DAD^-1,
1071:    where D is an appropriate diagonal matrix. This improves the accuracy of
1072:    the computed results in some cases. See the SLEPc Users Manual for details.

1074:    Balancing makes sense only for non-Hermitian problems when the required
1075:    precision is high (i.e. a small tolerance such as 1e-15).

1077:    By default, balancing is disabled. The two-sided method is much more
1078:    effective than the one-sided counterpart, but it requires the system
1079:    matrices to have the MatMultTranspose operation defined.

1081:    The parameter 'its' is the number of iterations performed by the method. The
1082:    cutoff value is used only in the two-side variant. Use PETSC_DEFAULT to assign
1083:    a reasonably good value.

1085:    User-defined balancing is allowed provided that the corresponding matrix
1086:    is set via STSetBalanceMatrix.

1088:    Level: intermediate

1090: .seealso: EPSGetBalance(), EPSBalance, STSetBalanceMatrix()
1091: @*/
1092: PetscErrorCode EPSSetBalance(EPS eps,EPSBalance bal,PetscInt its,PetscReal cutoff)
1093: {
1099:   switch (bal) {
1100:     case EPS_BALANCE_NONE:
1101:     case EPS_BALANCE_ONESIDE:
1102:     case EPS_BALANCE_TWOSIDE:
1103:     case EPS_BALANCE_USER:
1104:       if (eps->balance != bal) {
1105:         eps->state = EPS_STATE_INITIAL;
1106:         eps->balance = bal;
1107:       }
1108:       break;
1109:     default:
1110:       SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid value of argument 'bal'");
1111:   }
1112:   if (its==PETSC_DECIDE || its==PETSC_DEFAULT) eps->balance_its = 5;
1113:   else {
1114:     if (its <= 0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be > 0");
1115:     eps->balance_its = its;
1116:   }
1117:   if (cutoff==PETSC_DECIDE || cutoff==PETSC_DEFAULT) eps->balance_cutoff = 1e-8;
1118:   else {
1119:     if (cutoff <= 0.0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of cutoff. Must be > 0");
1120:     eps->balance_cutoff = cutoff;
1121:   }
1122:   return(0);
1123: }

1125: /*@C
1126:    EPSGetBalance - Gets the balancing type used by the EPS object, and the
1127:    associated parameters.

1129:    Not Collective

1131:    Input Parameter:
1132: .  eps - the eigensolver context

1134:    Output Parameters:
1135: +  bal    - the balancing method
1136: .  its    - number of iterations of the balancing algorithm
1137: -  cutoff - cutoff value

1139:    Level: intermediate

1141:    Note:
1142:    The user can specify NULL for any parameter that is not needed.

1144: .seealso: EPSSetBalance(), EPSBalance
1145: @*/
1146: PetscErrorCode EPSGetBalance(EPS eps,EPSBalance *bal,PetscInt *its,PetscReal *cutoff)
1147: {
1150:   if (bal)    *bal = eps->balance;
1151:   if (its)    *its = eps->balance_its;
1152:   if (cutoff) *cutoff = eps->balance_cutoff;
1153:   return(0);
1154: }

1156: /*@
1157:    EPSSetTwoSided - Sets the solver to use a two-sided variant so that left
1158:    eigenvectors are also computed.

1160:    Logically Collective on eps

1162:    Input Parameters:
1163: +  eps      - the eigensolver context
1164: -  twosided - whether the two-sided variant is to be used or not

1166:    Options Database Keys:
1167: .  -eps_two_sided <boolean> - Sets/resets the twosided flag

1169:    Notes:
1170:    If the user sets twosided=PETSC_TRUE then the solver uses a variant of
1171:    the algorithm that computes both right and left eigenvectors. This is
1172:    usually much more costly. This option is not available in all solvers.

1174:    When using two-sided solvers, the problem matrices must have both the
1175:    MatMult and MatMultTranspose operations defined.

1177:    Level: advanced

1179: .seealso: EPSGetTwoSided(), EPSGetLeftEigenvector()
1180: @*/
1181: PetscErrorCode EPSSetTwoSided(EPS eps,PetscBool twosided)
1182: {
1186:   if (twosided!=eps->twosided) {
1187:     eps->twosided = twosided;
1188:     eps->state    = EPS_STATE_INITIAL;
1189:   }
1190:   return(0);
1191: }

1193: /*@
1194:    EPSGetTwoSided - Returns the flag indicating whether a two-sided variant
1195:    of the algorithm is being used or not.

1197:    Not Collective

1199:    Input Parameter:
1200: .  eps - the eigensolver context

1202:    Output Parameter:
1203: .  twosided - the returned flag

1205:    Level: advanced

1207: .seealso: EPSSetTwoSided()
1208: @*/
1209: PetscErrorCode EPSGetTwoSided(EPS eps,PetscBool *twosided)
1210: {
1214:   *twosided = eps->twosided;
1215:   return(0);
1216: }

1218: /*@
1219:    EPSSetTrueResidual - Specifies if the solver must compute the true residual
1220:    explicitly or not.

1222:    Logically Collective on eps

1224:    Input Parameters:
1225: +  eps     - the eigensolver context
1226: -  trueres - whether true residuals are required or not

1228:    Options Database Keys:
1229: .  -eps_true_residual <boolean> - Sets/resets the boolean flag 'trueres'

1231:    Notes:
1232:    If the user sets trueres=PETSC_TRUE then the solver explicitly computes
1233:    the true residual for each eigenpair approximation, and uses it for
1234:    convergence testing. Computing the residual is usually an expensive
1235:    operation. Some solvers (e.g., Krylov solvers) can avoid this computation
1236:    by using a cheap estimate of the residual norm, but this may sometimes
1237:    give inaccurate results (especially if a spectral transform is being
1238:    used). On the contrary, preconditioned eigensolvers (e.g., Davidson solvers)
1239:    do rely on computing the true residual, so this option is irrelevant for them.

1241:    Level: advanced

1243: .seealso: EPSGetTrueResidual()
1244: @*/
1245: PetscErrorCode EPSSetTrueResidual(EPS eps,PetscBool trueres)
1246: {
1250:   eps->trueres = trueres;
1251:   return(0);
1252: }

1254: /*@
1255:    EPSGetTrueResidual - Returns the flag indicating whether true
1256:    residuals must be computed explicitly or not.

1258:    Not Collective

1260:    Input Parameter:
1261: .  eps - the eigensolver context

1263:    Output Parameter:
1264: .  trueres - the returned flag

1266:    Level: advanced

1268: .seealso: EPSSetTrueResidual()
1269: @*/
1270: PetscErrorCode EPSGetTrueResidual(EPS eps,PetscBool *trueres)
1271: {
1275:   *trueres = eps->trueres;
1276:   return(0);
1277: }

1279: /*@
1280:    EPSSetTrackAll - Specifies if the solver must compute the residual norm of all
1281:    approximate eigenpairs or not.

1283:    Logically Collective on eps

1285:    Input Parameters:
1286: +  eps      - the eigensolver context
1287: -  trackall - whether to compute all residuals or not

1289:    Notes:
1290:    If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
1291:    the residual norm for each eigenpair approximation. Computing the residual is
1292:    usually an expensive operation and solvers commonly compute only the residual
1293:    associated to the first unconverged eigenpair.

1295:    The option '-eps_monitor_all' automatically activates this option.

1297:    Level: developer

1299: .seealso: EPSGetTrackAll()
1300: @*/
1301: PetscErrorCode EPSSetTrackAll(EPS eps,PetscBool trackall)
1302: {
1306:   eps->trackall = trackall;
1307:   return(0);
1308: }

1310: /*@
1311:    EPSGetTrackAll - Returns the flag indicating whether all residual norms must
1312:    be computed or not.

1314:    Not Collective

1316:    Input Parameter:
1317: .  eps - the eigensolver context

1319:    Output Parameter:
1320: .  trackall - the returned flag

1322:    Level: developer

1324: .seealso: EPSSetTrackAll()
1325: @*/
1326: PetscErrorCode EPSGetTrackAll(EPS eps,PetscBool *trackall)
1327: {
1331:   *trackall = eps->trackall;
1332:   return(0);
1333: }

1335: /*@
1336:    EPSSetPurify - Deactivate eigenvector purification (which is activated by default).

1338:    Logically Collective on eps

1340:    Input Parameters:
1341: +  eps    - the eigensolver context
1342: -  purify - whether purification is required or not

1344:    Options Database Keys:
1345: .  -eps_purify <boolean> - Sets/resets the boolean flag 'purify'

1347:    Notes:
1348:    By default, eigenvectors of generalized symmetric eigenproblems are purified
1349:    in order to purge directions in the nullspace of matrix B. If the user knows
1350:    that B is non-singular, then purification can be safely deactivated and some
1351:    computational cost is avoided (this is particularly important in interval computations).

1353:    Level: intermediate

1355: .seealso: EPSGetPurify(), EPSSetInterval()
1356: @*/
1357: PetscErrorCode EPSSetPurify(EPS eps,PetscBool purify)
1358: {
1362:   if (purify!=eps->purify) {
1363:     eps->purify = purify;
1364:     eps->state  = EPS_STATE_INITIAL;
1365:   }
1366:   return(0);
1367: }

1369: /*@
1370:    EPSGetPurify - Returns the flag indicating whether purification is activated
1371:    or not.

1373:    Not Collective

1375:    Input Parameter:
1376: .  eps - the eigensolver context

1378:    Output Parameter:
1379: .  purify - the returned flag

1381:    Level: intermediate

1383: .seealso: EPSSetPurify()
1384: @*/
1385: PetscErrorCode EPSGetPurify(EPS eps,PetscBool *purify)
1386: {
1390:   *purify = eps->purify;
1391:   return(0);
1392: }

1394: /*@C
1395:    EPSSetOptionsPrefix - Sets the prefix used for searching for all
1396:    EPS options in the database.

1398:    Logically Collective on eps

1400:    Input Parameters:
1401: +  eps - the eigensolver context
1402: -  prefix - the prefix string to prepend to all EPS option requests

1404:    Notes:
1405:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1406:    The first character of all runtime options is AUTOMATICALLY the
1407:    hyphen.

1409:    For example, to distinguish between the runtime options for two
1410:    different EPS contexts, one could call
1411: .vb
1412:       EPSSetOptionsPrefix(eps1,"eig1_")
1413:       EPSSetOptionsPrefix(eps2,"eig2_")
1414: .ve

1416:    Level: advanced

1418: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
1419: @*/
1420: PetscErrorCode EPSSetOptionsPrefix(EPS eps,const char *prefix)
1421: {

1426:   if (!eps->st) { EPSGetST(eps,&eps->st); }
1427:   STSetOptionsPrefix(eps->st,prefix);
1428:   if (!eps->V) { EPSGetBV(eps,&eps->V); }
1429:   BVSetOptionsPrefix(eps->V,prefix);
1430:   if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
1431:   DSSetOptionsPrefix(eps->ds,prefix);
1432:   if (!eps->rg) { EPSGetRG(eps,&eps->rg); }
1433:   RGSetOptionsPrefix(eps->rg,prefix);
1434:   PetscObjectSetOptionsPrefix((PetscObject)eps,prefix);
1435:   return(0);
1436: }

1438: /*@C
1439:    EPSAppendOptionsPrefix - Appends to the prefix used for searching for all
1440:    EPS options in the database.

1442:    Logically Collective on eps

1444:    Input Parameters:
1445: +  eps - the eigensolver context
1446: -  prefix - the prefix string to prepend to all EPS option requests

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

1452:    Level: advanced

1454: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
1455: @*/
1456: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,const char *prefix)
1457: {

1462:   if (!eps->st) { EPSGetST(eps,&eps->st); }
1463:   STAppendOptionsPrefix(eps->st,prefix);
1464:   if (!eps->V) { EPSGetBV(eps,&eps->V); }
1465:   BVAppendOptionsPrefix(eps->V,prefix);
1466:   if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
1467:   DSAppendOptionsPrefix(eps->ds,prefix);
1468:   if (!eps->rg) { EPSGetRG(eps,&eps->rg); }
1469:   RGAppendOptionsPrefix(eps->rg,prefix);
1470:   PetscObjectAppendOptionsPrefix((PetscObject)eps,prefix);
1471:   return(0);
1472: }

1474: /*@C
1475:    EPSGetOptionsPrefix - Gets the prefix used for searching for all
1476:    EPS options in the database.

1478:    Not Collective

1480:    Input Parameters:
1481: .  eps - the eigensolver context

1483:    Output Parameters:
1484: .  prefix - pointer to the prefix string used is returned

1486:    Note:
1487:    On the Fortran side, the user should pass in a string 'prefix' of
1488:    sufficient length to hold the prefix.

1490:    Level: advanced

1492: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
1493: @*/
1494: PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])
1495: {

1501:   PetscObjectGetOptionsPrefix((PetscObject)eps,prefix);
1502:   return(0);
1503: }