Actual source code: opts.c

  1: /*
  2:       EPS 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-2011, Universitat Politecnica de Valencia, Spain

  9:    This file is part of SLEPc.
 10:       
 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 <private/epsimpl.h>   /*I "slepceps.h" I*/

 29: /*@
 30:    EPSSetFromOptions - Sets EPS options from the options database.
 31:    This routine must be called before EPSSetUp() if the user is to be 
 32:    allowed to set the solver type. 

 34:    Collective on EPS

 36:    Input Parameters:
 37: .  eps - the eigensolver context

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

 42:    Level: beginner
 43: @*/
 44: PetscErrorCode EPSSetFromOptions(EPS eps)
 45: {
 46:   PetscErrorCode   ierr;
 47:   char             type[256],monfilename[PETSC_MAX_PATH_LEN];
 48:   PetscBool        flg,val;
 49:   PetscReal        r,nrma,nrmb,array[2]={0,0};
 50:   PetscScalar      s;
 51:   PetscInt         i,j,k;
 52:   const char       *bal_list[4] = {"none","oneside","twoside","user"};
 53:   PetscViewer      monviewer;
 54:   SlepcConvMonitor ctx;

 58:   if (!EPSRegisterAllCalled) { EPSRegisterAll(PETSC_NULL); }
 59:   PetscOptionsBegin(((PetscObject)eps)->comm,((PetscObject)eps)->prefix,"Eigenproblem Solver (EPS) Options","EPS");
 60:     PetscOptionsList("-eps_type","Eigenproblem Solver method","EPSSetType",EPSList,(char*)(((PetscObject)eps)->type_name?((PetscObject)eps)->type_name:EPSKRYLOVSCHUR),type,256,&flg);
 61:     if (flg) {
 62:       EPSSetType(eps,type);
 63:     }
 64:     /*
 65:       Set the type if it was never set.
 66:     */
 67:     if (!((PetscObject)eps)->type_name) {
 68:       EPSSetType(eps,EPSKRYLOVSCHUR);
 69:     }

 71:     PetscOptionsBoolGroupBegin("-eps_hermitian","hermitian eigenvalue problem","EPSSetProblemType",&flg);
 72:     if (flg) {EPSSetProblemType(eps,EPS_HEP);}
 73:     PetscOptionsBoolGroup("-eps_gen_hermitian","generalized hermitian eigenvalue problem","EPSSetProblemType",&flg);
 74:     if (flg) {EPSSetProblemType(eps,EPS_GHEP);}
 75:     PetscOptionsBoolGroup("-eps_non_hermitian","non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
 76:     if (flg) {EPSSetProblemType(eps,EPS_NHEP);}
 77:     PetscOptionsBoolGroup("-eps_gen_non_hermitian","generalized non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
 78:     if (flg) {EPSSetProblemType(eps,EPS_GNHEP);}
 79:     PetscOptionsBoolGroupEnd("-eps_pos_gen_non_hermitian","generalized non-hermitian eigenvalue problem with positive semi-definite B","EPSSetProblemType",&flg);
 80:     if (flg) {EPSSetProblemType(eps,EPS_PGNHEP);}

 82:     PetscOptionsBoolGroupBegin("-eps_ritz","Rayleigh-Ritz extraction","EPSSetExtraction",&flg);
 83:     if (flg) {EPSSetExtraction(eps,EPS_RITZ);}
 84:     PetscOptionsBoolGroup("-eps_harmonic","harmonic Ritz extraction","EPSSetExtraction",&flg);
 85:     if (flg) {EPSSetExtraction(eps,EPS_HARMONIC);}
 86:     PetscOptionsBoolGroup("-eps_harmonic_relative","relative harmonic Ritz extraction","EPSSetExtraction",&flg);
 87:     if (flg) {EPSSetExtraction(eps,EPS_HARMONIC_RELATIVE);}
 88:     PetscOptionsBoolGroup("-eps_harmonic_right","right harmonic Ritz extraction","EPSSetExtraction",&flg);
 89:     if (flg) {EPSSetExtraction(eps,EPS_HARMONIC_RIGHT);}
 90:     PetscOptionsBoolGroup("-eps_harmonic_largest","largest harmonic Ritz extraction","EPSSetExtraction",&flg);
 91:     if (flg) {EPSSetExtraction(eps,EPS_HARMONIC_LARGEST);}
 92:     PetscOptionsBoolGroup("-eps_refined","refined Ritz extraction","EPSSetExtraction",&flg);
 93:     if (flg) {EPSSetExtraction(eps,EPS_REFINED);}
 94:     PetscOptionsBoolGroupEnd("-eps_refined_harmonic","refined harmonic Ritz extraction","EPSSetExtraction",&flg);
 95:     if (flg) {EPSSetExtraction(eps,EPS_REFINED_HARMONIC);}

 97:     if (!eps->balance) eps->balance = EPS_BALANCE_NONE;
 98:     PetscOptionsEList("-eps_balance","Balancing method","EPSSetBalance",bal_list,4,bal_list[eps->balance-EPS_BALANCE_NONE],&i,&flg);
 99:     if (flg) { eps->balance = (EPSBalance)(i+EPS_BALANCE_NONE); }
100:     r = j = PETSC_IGNORE;
101:     PetscOptionsInt("-eps_balance_its","Number of iterations in balancing","EPSSetBalance",eps->balance_its,&j,PETSC_NULL);
102:     PetscOptionsReal("-eps_balance_cutoff","Cutoff value in balancing","EPSSetBalance",eps->balance_cutoff,&r,PETSC_NULL);
103:     EPSSetBalance(eps,(EPSBalance)PETSC_IGNORE,j,r);

105:     r = i = PETSC_IGNORE;
106:     PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&i,PETSC_NULL);
107:     PetscOptionsReal("-eps_tol","Tolerance","EPSSetTolerances",eps->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:eps->tol,&r,PETSC_NULL);
108:     EPSSetTolerances(eps,r,i);
109:     PetscOptionsBoolGroupBegin("-eps_conv_eig","relative error convergence test","EPSSetConvergenceTest",&flg);
110:     if (flg) {EPSSetConvergenceTest(eps,EPS_CONV_EIG);}
111:     PetscOptionsBoolGroup("-eps_conv_norm","Convergence test relative to the eigenvalue and the matrix norms","EPSSetConvergenceTest",&flg);
112:     if (flg) {EPSSetConvergenceTest(eps,EPS_CONV_NORM);}
113:     PetscOptionsBoolGroupEnd("-eps_conv_abs","Absolute error convergence test","EPSSetConvergenceTest",&flg);
114:     if (flg) {EPSSetConvergenceTest(eps,EPS_CONV_ABS);}

116:     i = j = k = PETSC_IGNORE;
117:     PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&i,PETSC_NULL);
118:     PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&j,PETSC_NULL);
119:     PetscOptionsInt("-eps_mpd","Maximum dimension of projected problem","EPSSetDimensions",eps->mpd,&k,PETSC_NULL);
120:     EPSSetDimensions(eps,i,j,k);
121: 
122:     /* -----------------------------------------------------------------------*/
123:     /*
124:       Cancels all monitors hardwired into code before call to EPSSetFromOptions()
125:     */
126:     flg  = PETSC_FALSE;
127:     PetscOptionsBool("-eps_monitor_cancel","Remove any hardwired monitor routines","EPSMonitorCancel",flg,&flg,PETSC_NULL);
128:     if (flg) {
129:       EPSMonitorCancel(eps);
130:     }
131:     /*
132:       Prints approximate eigenvalues and error estimates at each iteration
133:     */
134:     PetscOptionsString("-eps_monitor","Monitor first unconverged approximate eigenvalue and error estimate","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
135:     if (flg) {
136:       PetscViewerASCIIOpen(((PetscObject)eps)->comm,monfilename,&monviewer);
137:       EPSMonitorSet(eps,EPSMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
138:     }
139:     PetscOptionsString("-eps_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
140:     if (flg) {
141:       PetscNew(struct _n_SlepcConvMonitor,&ctx);
142:       PetscViewerASCIIOpen(((PetscObject)eps)->comm,monfilename,&ctx->viewer);
143:       EPSMonitorSet(eps,EPSMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
144:     }
145:     PetscOptionsString("-eps_monitor_all","Monitor approximate eigenvalues and error estimates","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
146:     if (flg) {
147:       PetscViewerASCIIOpen(((PetscObject)eps)->comm,monfilename,&monviewer);
148:       EPSMonitorSet(eps,EPSMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
149:       EPSSetTrackAll(eps,PETSC_TRUE);
150:     }
151:     flg = PETSC_FALSE;
152:     PetscOptionsBool("-eps_monitor_draw","Monitor first unconverged approximate eigenvalue and error estimate graphically","EPSMonitorSet",flg,&flg,PETSC_NULL);
153:     if (flg) {
154:       EPSMonitorSet(eps,EPSMonitorLG,PETSC_NULL,PETSC_NULL);
155:     }
156:     flg = PETSC_FALSE;
157:     PetscOptionsBool("-eps_monitor_draw_all","Monitor error estimates graphically","EPSMonitorSet",flg,&flg,PETSC_NULL);
158:     if (flg) {
159:       EPSMonitorSet(eps,EPSMonitorLGAll,PETSC_NULL,PETSC_NULL);
160:       EPSSetTrackAll(eps,PETSC_TRUE);
161:     }
162:   /* -----------------------------------------------------------------------*/
163:     PetscOptionsScalar("-eps_target","Value of the target","EPSSetTarget",eps->target,&s,&flg);
164:     if (flg) {
165:       EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);
166:       EPSSetTarget(eps,s);
167:     }
168:     k = 2;
169:     PetscOptionsRealArray("-eps_interval","Computational interval (two real values separated with a comma without spaces)","EPSSetInterval",array,&k,&flg);
170:     if (flg) {
171:       if (k<2) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_SIZ,"Must pass two values in -eps_interval (comma-separated without spaces)");
172:       EPSSetWhichEigenpairs(eps,EPS_ALL);
173:       EPSSetInterval(eps,array[0],array[1]);
174:     }

176:     PetscOptionsBoolGroupBegin("-eps_largest_magnitude","compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
177:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE);}
178:     PetscOptionsBoolGroup("-eps_smallest_magnitude","compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
179:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE);}
180:     PetscOptionsBoolGroup("-eps_largest_real","compute largest real parts","EPSSetWhichEigenpairs",&flg);
181:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL);}
182:     PetscOptionsBoolGroup("-eps_smallest_real","compute smallest real parts","EPSSetWhichEigenpairs",&flg);
183:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);}
184:     PetscOptionsBoolGroup("-eps_largest_imaginary","compute largest imaginary parts","EPSSetWhichEigenpairs",&flg);
185:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY);}
186:     PetscOptionsBoolGroup("-eps_smallest_imaginary","compute smallest imaginary parts","EPSSetWhichEigenpairs",&flg);
187:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY);}
188:     PetscOptionsBoolGroup("-eps_target_magnitude","compute nearest eigenvalues to target","EPSSetWhichEigenpairs",&flg);
189:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);}
190:     PetscOptionsBoolGroup("-eps_target_real","compute eigenvalues with real parts close to target","EPSSetWhichEigenpairs",&flg);
191:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_TARGET_REAL);}
192:     PetscOptionsBoolGroup("-eps_target_imaginary","compute eigenvalues with imaginary parts close to target","EPSSetWhichEigenpairs",&flg);
193:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_TARGET_IMAGINARY);}
194:     PetscOptionsBoolGroupEnd("-eps_all","compute all eigenvalues in an interval","EPSSetWhichEigenpairs",&flg);
195:     if (flg) {EPSSetWhichEigenpairs(eps,EPS_ALL);}

197:     PetscOptionsBool("-eps_left_vectors","Compute left eigenvectors also","EPSSetLeftVectorsWanted",eps->leftvecs,&val,&flg);
198:     if (flg) {
199:       EPSSetLeftVectorsWanted(eps,val);
200:     }
201:     PetscOptionsBool("-eps_true_residual","Compute true residuals explicitly","EPSSetTrueResidual",eps->trueres,&val,&flg);
202:     if (flg) {
203:       EPSSetTrueResidual(eps,val);
204:     }

206:     nrma = nrmb = PETSC_IGNORE;
207:     PetscOptionsReal("-eps_norm_a","Norm of matrix A","EPSSetMatrixNorms",eps->nrma,&nrma,PETSC_NULL);
208:     PetscOptionsReal("-eps_norm_b","Norm of matrix B","EPSSetMatrixNorms",eps->nrmb,&nrmb,PETSC_NULL);
209:     EPSSetMatrixNorms(eps,nrma,nrmb,eps->adaptive);
210:     PetscOptionsBool("-eps_norms_adaptive","Update the value of matrix norms adaptively","EPSSetMatrixNorms",eps->adaptive,&val,&flg);
211:     if (flg) {
212:       EPSSetMatrixNorms(eps,PETSC_IGNORE,PETSC_IGNORE,val);
213:     }

215:     PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",0);
216:     PetscOptionsName("-eps_view_binary","Save the matrices associated to the eigenproblem","EPSSetFromOptions",0);
217:     PetscOptionsName("-eps_plot_eigs","Make a plot of the computed eigenvalues","EPSSolve",0);
218: 
219:     if (eps->ops->setfromoptions) {
220:       (*eps->ops->setfromoptions)(eps);
221:     }
222:     PetscObjectProcessOptionsHandlers((PetscObject)eps);
223:   PetscOptionsEnd();

225:   if (!eps->ip) { EPSGetIP(eps,&eps->ip); }
226:   IPSetFromOptions(eps->ip);
227:   STSetFromOptions(eps->OP);
228:   PetscRandomSetFromOptions(eps->rand);
229:   return(0);
230: }

234: /*@
235:    EPSGetTolerances - Gets the tolerance and maximum iteration count used
236:    by the EPS convergence tests. 

238:    Not Collective

240:    Input Parameter:
241: .  eps - the eigensolver context
242:   
243:    Output Parameters:
244: +  tol - the convergence tolerance
245: -  maxits - maximum number of iterations

247:    Notes:
248:    The user can specify PETSC_NULL for any parameter that is not needed.

250:    Level: intermediate

252: .seealso: EPSSetTolerances()
253: @*/
254: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,PetscInt *maxits)
255: {
258:   if (tol)    *tol    = eps->tol;
259:   if (maxits) *maxits = eps->max_it;
260:   return(0);
261: }

265: /*@
266:    EPSSetTolerances - Sets the tolerance and maximum iteration count used
267:    by the EPS convergence tests. 

269:    Logically Collective on EPS

271:    Input Parameters:
272: +  eps - the eigensolver context
273: .  tol - the convergence tolerance
274: -  maxits - maximum number of iterations to use

276:    Options Database Keys:
277: +  -eps_tol <tol> - Sets the convergence tolerance
278: -  -eps_max_it <maxits> - Sets the maximum number of iterations allowed

280:    Notes:
281:    Use PETSC_IGNORE for an argument that need not be changed.

283:    Use PETSC_DECIDE for maxits to assign a reasonably good value, which is 
284:    dependent on the solution method.

286:    Level: intermediate

288: .seealso: EPSGetTolerances()
289: @*/
290: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,PetscInt maxits)
291: {
296:   if (tol != PETSC_IGNORE) {
297:     if (tol == PETSC_DEFAULT) {
298:       eps->tol = PETSC_DEFAULT;
299:     } else {
300:       if (tol < 0.0) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
301:       eps->tol = tol;
302:     }
303:   }
304:   if (maxits != PETSC_IGNORE) {
305:     if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
306:       eps->max_it = 0;
307:       eps->setupcalled = 0;
308:     } else {
309:       if (maxits < 0) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
310:       eps->max_it = maxits;
311:     }
312:   }
313:   return(0);
314: }

318: /*@
319:    EPSGetDimensions - Gets the number of eigenvalues to compute
320:    and the dimension of the subspace.

322:    Not Collective

324:    Input Parameter:
325: .  eps - the eigensolver context
326:   
327:    Output Parameters:
328: +  nev - number of eigenvalues to compute
329: .  ncv - the maximum dimension of the subspace to be used by the solver
330: -  mpd - the maximum dimension allowed for the projected problem

332:    Notes:
333:    The user can specify PETSC_NULL for any parameter that is not needed.

335:    Level: intermediate

337: .seealso: EPSSetDimensions()
338: @*/
339: PetscErrorCode EPSGetDimensions(EPS eps,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
340: {
343:   if (nev) *nev = eps->nev;
344:   if (ncv) *ncv = eps->ncv;
345:   if (mpd) *mpd = eps->mpd;
346:   return(0);
347: }

351: /*@
352:    EPSSetDimensions - Sets the number of eigenvalues to compute
353:    and the dimension of the subspace.

355:    Logically Collective on EPS

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

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

368:    Notes:
369:    Use PETSC_IGNORE to retain the previous value of any parameter.

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

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

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

383:    When computing all eigenvalues in an interval, see EPSSetInterval(), the
384:    meaning of nev changes. In that case, the number of eigenvalues in the
385:    interval is not known a priori; the meaning of nev is then the number of
386:    eigenvalues that are computed at a time when sweeping the interval from one
387:    end to the other. The value of nev in this case may have an impact on overall
388:    performance. The value of ncv should not be assigned in this case.

390:    Level: intermediate

392: .seealso: EPSGetDimensions(), EPSSetInterval()
393: @*/
394: PetscErrorCode EPSSetDimensions(EPS eps,PetscInt nev,PetscInt ncv,PetscInt mpd)
395: {
401:   if (nev != PETSC_IGNORE) {
402:     if (nev<1) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
403:     eps->nev = nev;
404:     eps->setupcalled = 0;
405:   }
406:   if (ncv != PETSC_IGNORE) {
407:     if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
408:       eps->ncv = 0;
409:     } else {
410:       if (ncv<1) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
411:       eps->ncv = ncv;
412:     }
413:     eps->setupcalled = 0;
414:   }
415:   if (mpd != PETSC_IGNORE) {
416:     if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
417:       eps->mpd = 0;
418:     } else {
419:       if (mpd<1) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
420:       eps->mpd = mpd;
421:     }
422:   }
423:   return(0);
424: }

428: /*@
429:    EPSSetWhichEigenpairs - Specifies which portion of the spectrum is 
430:    to be sought.

432:    Logically Collective on EPS

434:    Input Parameters:
435: +  eps   - eigensolver context obtained from EPSCreate()
436: -  which - the portion of the spectrum to be sought

438:    Possible values:
439:    The parameter 'which' can have one of these values
440:     
441: +     EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
442: .     EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
443: .     EPS_LARGEST_REAL - largest real parts
444: .     EPS_SMALLEST_REAL - smallest real parts
445: .     EPS_LARGEST_IMAGINARY - largest imaginary parts
446: .     EPS_SMALLEST_IMAGINARY - smallest imaginary parts
447: .     EPS_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
448: .     EPS_TARGET_REAL - eigenvalues with real part closest to target
449: .     EPS_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
450: .     EPS_ALL - all eigenvalues contained in a given interval
451: -     EPS_WHICH_USER - user defined ordering set with EPSSetEigenvalueComparison()

453:    Options Database Keys:
454: +   -eps_largest_magnitude - Sets largest eigenvalues in magnitude
455: .   -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
456: .   -eps_largest_real - Sets largest real parts
457: .   -eps_smallest_real - Sets smallest real parts
458: .   -eps_largest_imaginary - Sets largest imaginary parts
459: .   -eps_smallest_imaginary - Sets smallest imaginary parts
460: .   -eps_target_magnitude - Sets eigenvalues closest to target
461: .   -eps_target_real - Sets real parts closest to target
462: .   -eps_target_imaginary - Sets imaginary parts closest to target
463: -   -eps_all - Sets all eigenvalues in an interval

465:    Notes:
466:    Not all eigensolvers implemented in EPS account for all the possible values
467:    stated above. Also, some values make sense only for certain types of 
468:    problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY
469:    and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part 
470:    for eigenvalue selection.
471:     
472:    The target is a scalar value provided with EPSSetTarget().

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

477:    EPS_ALL is intended for use in combination with an interval (see
478:    EPSSetInterval()), when all eigenvalues within the interval are requested.
479:    In that case, the number of eigenvalues is unknown, so the nev parameter
480:    has a different sense, see EPSSetDimensions().

482:    Level: intermediate

484: .seealso: EPSGetWhichEigenpairs(), EPSSetTarget(), EPSSetInterval(),
485:           EPSSetDimensions(), EPSSetEigenvalueComparison(),
486:           EPSSortEigenvalues(), EPSWhich
487: @*/
488: PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)
489: {
493:   if (which!=PETSC_IGNORE) {
494:     if (which==PETSC_DECIDE || which==PETSC_DEFAULT) eps->which = (EPSWhich)0;
495:     else switch (which) {
496:       case EPS_LARGEST_MAGNITUDE:
497:       case EPS_SMALLEST_MAGNITUDE:
498:       case EPS_LARGEST_REAL:
499:       case EPS_SMALLEST_REAL:
500:       case EPS_LARGEST_IMAGINARY:
501:       case EPS_SMALLEST_IMAGINARY:
502:       case EPS_TARGET_MAGNITUDE:
503:       case EPS_TARGET_REAL:
504: #if defined(PETSC_USE_COMPLEX)
505:       case EPS_TARGET_IMAGINARY:
506: #endif
507:       case EPS_ALL:
508:       case EPS_WHICH_USER:
509:         if (eps->which != which) {
510:           eps->setupcalled = 0;
511:           eps->which = which;
512:         }
513:         break;
514:       default:
515:         SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
516:     }
517:   }
518:   return(0);
519: }

523: /*@C
524:    EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be 
525:    sought.

527:    Not Collective

529:    Input Parameter:
530: .  eps - eigensolver context obtained from EPSCreate()

532:    Output Parameter:
533: .  which - the portion of the spectrum to be sought

535:    Notes:
536:    See EPSSetWhichEigenpairs() for possible values of 'which'.

538:    Level: intermediate

540: .seealso: EPSSetWhichEigenpairs(), EPSWhich
541: @*/
542: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)
543: {
547:   *which = eps->which;
548:   return(0);
549: }

553: /*@
554:    EPSSetLeftVectorsWanted - Specifies which eigenvectors are required.

556:    Logically Collective on EPS

558:    Input Parameters:
559: +  eps      - the eigensolver context
560: -  leftvecs - whether left eigenvectors are required or not

562:    Options Database Keys:
563: .  -eps_left_vectors <boolean> - Sets/resets the boolean flag 'leftvecs'

565:    Notes:
566:    If the user sets leftvecs=PETSC_TRUE then the solver uses a variant of
567:    the algorithm that computes both right and left eigenvectors. This is
568:    usually much more costly. This option is not available in all solvers.

570:    Level: intermediate

572: .seealso: EPSGetLeftVectorsWanted(), EPSGetEigenvectorLeft()
573: @*/
574: PetscErrorCode EPSSetLeftVectorsWanted(EPS eps,PetscBool leftvecs)
575: {
579:   if (eps->leftvecs != leftvecs) {
580:     eps->leftvecs = leftvecs;
581:     eps->setupcalled = 0;
582:   }
583:   return(0);
584: }

588: /*@
589:    EPSGetLeftVectorsWanted - Returns the flag indicating whether left 
590:    eigenvectors are required or not.

592:    Not Collective

594:    Input Parameter:
595: .  eps - the eigensolver context

597:    Output Parameter:
598: .  leftvecs - the returned flag

600:    Level: intermediate

602: .seealso: EPSSetLeftVectorsWanted(), EPSWhich
603: @*/
604: PetscErrorCode EPSGetLeftVectorsWanted(EPS eps,PetscBool *leftvecs)
605: {
609:   *leftvecs = eps->leftvecs;
610:   return(0);
611: }

615: /*@
616:    EPSSetMatrixNorms - Gives the reference values of the matrix norms
617:    and specifies whether these values should be improved adaptively.

619:    Logically Collective on EPS

621:    Input Parameters:
622: +  eps      - the eigensolver context
623: .  nrma     - a reference value for the norm of matrix A
624: .  nrmb     - a reference value for the norm of matrix B
625: -  adaptive - whether matrix norms are improved adaptively

627:    Options Database Keys:
628: +  -eps_norm_a <nrma> - norm of A
629: .  -eps_norm_b <nrma> - norm of B
630: -  -eps_norms_adaptive <boolean> - Sets/resets the boolean flag 'adaptive'

632:    Notes:
633:    If the user sets adaptive=PETSC_FALSE then the solver uses the values
634:    of nrma and nrmb for the matrix norms, and these values do not change
635:    throughout the iteration.

637:    If the user sets adaptive=PETSC_TRUE then the solver tries to adaptively
638:    improve the supplied values, with the numerical information generated
639:    during the iteration. This option is not available in all solvers.

641:    If a passed value is PETSC_DEFAULT, the corresponding norm will be set to 1.
642:    If a passed value is PETSC_DETERMINE, the corresponding norm will be computed
643:    as the NORM_INFINITY with MatNorm().

645:    Level: intermediate

647: .seealso: EPSGetMatrixNorms()
648: @*/
649: PetscErrorCode EPSSetMatrixNorms(EPS eps,PetscReal nrma,PetscReal nrmb,PetscBool adaptive)
650: {
656:   if (nrma != PETSC_IGNORE) {
657:     if (nrma == PETSC_DEFAULT) eps->nrma = 1.0;
658:     else if (nrma == PETSC_DETERMINE) {
659:       eps->nrma = nrma;
660:       eps->setupcalled = 0;
661:     } else {
662:       if (nrma < 0.0) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nrma. Must be > 0");
663:       eps->nrma = nrma;
664:     }
665:   }
666:   if (nrmb != PETSC_IGNORE) {
667:     if (!eps->isgeneralized) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_WRONG,"Norm of B only allowed in generalized problems");
668:     if (nrmb == PETSC_DEFAULT) eps->nrmb = 1.0;
669:     else if (nrmb == PETSC_DETERMINE) {
670:       eps->nrmb = nrmb;
671:       eps->setupcalled = 0;
672:     } else {
673:       if (nrmb < 0.0) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nrmb. Must be > 0");
674:       eps->nrmb = nrmb;
675:     }
676:   }
677:   if (eps->adaptive != adaptive) {
678:     eps->adaptive = adaptive;
679:     eps->setupcalled = 0;
680:     if (adaptive) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_SUP,"Sorry, adaptive norms are not implemented in this release.");
681:   }
682:   return(0);
683: }

687: /*@
688:    EPSGetMatrixNorms - Returns the value of the matrix norms (either set
689:    by the user or estimated by the solver) and the flag indicating whether
690:    the norms are being adaptively improved.

692:    Not Collective

694:    Input Parameter:
695: .  eps - the eigensolver context

697:    Output Parameters:
698: +  nrma     - the norm of matrix A
699: .  nrmb     - the norm of matrix B
700: -  adaptive - whether matrix norms are improved adaptively

702:    Level: intermediate

704: .seealso: EPSSetMatrixNorms()
705: @*/
706: PetscErrorCode EPSGetMatrixNorms(EPS eps,PetscReal *nrma,PetscReal *nrmb,PetscBool *adaptive)
707: {
710:   if (nrma) *nrma = eps->nrma;
711:   if (nrmb) *nrmb = eps->nrmb;
712:   if (adaptive) *adaptive = eps->adaptive;
713:   return(0);
714: }

718: /*@C
719:    EPSSetEigenvalueComparison - Specifies the eigenvalue comparison function
720:    when EPSSetWhichEigenpairs() is set to EPS_WHICH_USER.

722:    Logically Collective on EPS

724:    Input Parameters:
725: +  eps  - eigensolver context obtained from EPSCreate()
726: .  func - a pointer to the comparison function
727: -  ctx  - a context pointer (the last parameter to the comparison function)

729:    Calling Sequence of func:
730: $   func(EPS eps,PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)

732: +   eps    - eigensolver context obtained from EPSCreate()
733: .   ar     - real part of the 1st eigenvalue
734: .   ai     - imaginary part of the 1st eigenvalue
735: .   br     - real part of the 2nd eigenvalue
736: .   bi     - imaginary part of the 2nd eigenvalue
737: .   res    - result of comparison
738: -   ctx    - optional context, as set by EPSSetEigenvalueComparison()

740:    Note:
741:    The returning parameter 'res' can be:
742: +  negative - if the 1st eigenvalue is preferred to the 2st one
743: .  zero     - if both eigenvalues are equally preferred
744: -  positive - if the 2st eigenvalue is preferred to the 1st one

746:    Level: advanced

748: .seealso: EPSSetWhichEigenpairs(), EPSSortEigenvalues(), EPSWhich
749: @*/
750: PetscErrorCode EPSSetEigenvalueComparison(EPS eps,PetscErrorCode (*func)(EPS,PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)
751: {
754:   eps->which_func = func;
755:   eps->which_ctx  = ctx;
756:   eps->which      = EPS_WHICH_USER;
757:   return(0);
758: }

762: /*@C
763:    EPSSetConvergenceTestFunction - Sets a function to compute the error estimate
764:    used in the convergence test.

766:    Logically Collective on EPS

768:    Input Parameters:
769: +  eps  - eigensolver context obtained from EPSCreate()
770: .  func - a pointer to the convergence test function
771: -  ctx  - a context pointer (the last parameter to the convergence test function)

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

776: +   eps    - eigensolver context obtained from EPSCreate()
777: .   eigr   - real part of the eigenvalue
778: .   eigi   - imaginary part of the eigenvalue
779: .   res    - residual norm associated to the eigenpair
780: .   errest - (output) computed error estimate
781: -   ctx    - optional context, as set by EPSSetConvergenceTest()

783:    Note:
784:    If the error estimate returned by the convergence test function is less than
785:    the tolerance, then the eigenvalue is accepted as converged.

787:    Level: advanced

789: .seealso: EPSSetConvergenceTest(),EPSSetTolerances()
790: @*/
792: {
795:   eps->conv_func = func;
796:   eps->conv_ctx = ctx;
797:   if (func == EPSEigRelativeConverged) eps->conv = EPS_CONV_EIG;
798:   else if (func == EPSNormRelativeConverged) eps->conv = EPS_CONV_NORM;
799:   else if (func == EPSAbsoluteConverged) eps->conv = EPS_CONV_ABS;
800:   else eps->conv = EPS_CONV_USER;
801:   return(0);
802: }

806: /*@
807:    EPSSetConvergenceTest - Specifies how to compute the error estimate
808:    used in the convergence test.

810:    Logically Collective on EPS

812:    Input Parameters:
813: +  eps   - eigensolver context obtained from EPSCreate()
814: -  conv  - the type of convergence test

816:    Options Database Keys:
817: +  -eps_conv_abs - Sets the absolute convergence test
818: .  -eps_conv_eig - Sets the convergence test relative to the eigenvalue
819: -  -eps_conv_norm - Sets the convergence test relative to the matrix norms

821:    Note:
822:    The parameter 'conv' can have one of these values
823: +     EPS_CONV_ABS - absolute error ||r||
824: .     EPS_CONV_EIG - error relative to the eigenvalue l, ||r||/|l|
825: .     EPS_CONV_NORM - error relative to the matrix norms, ||r||/(||A||+|l|*||B||)
826: -     EPS_CONV_USER - function set by EPSSetConvergenceTestFunction() 

828:    Level: intermediate

830: .seealso: EPSGetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSConv
831: @*/
832: PetscErrorCode EPSSetConvergenceTest(EPS eps,EPSConv conv)
833: {
837:   switch(conv) {
838:   case EPS_CONV_EIG: eps->conv_func = EPSEigRelativeConverged; break;
839:   case EPS_CONV_NORM: eps->conv_func = EPSNormRelativeConverged; break;
840:   case EPS_CONV_ABS: eps->conv_func = EPSAbsoluteConverged; break;
841:   case EPS_CONV_USER: break;
842:   default:
843:     SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
844:   }
845:   eps->conv = conv;
846:   return(0);
847: }

851: /*@
852:    EPSGetConvergenceTest - Gets the method used to compute the error estimate
853:    used in the convergence test.

855:    Not Collective

857:    Input Parameters:
858: .  eps   - eigensolver context obtained from EPSCreate()

860:    Output Parameters:
861: .  conv  - the type of convergence test

863:    Level: intermediate

865: .seealso: EPSSetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSConv
866: @*/
867: PetscErrorCode EPSGetConvergenceTest(EPS eps,EPSConv *conv)
868: {
872:   *conv = eps->conv;
873:   return(0);
874: }

878: /*@
879:    EPSSetProblemType - Specifies the type of the eigenvalue problem.

881:    Logically Collective on EPS

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

887:    Options Database Keys:
888: +  -eps_hermitian - Hermitian eigenvalue problem
889: .  -eps_gen_hermitian - generalized Hermitian eigenvalue problem
890: .  -eps_non_hermitian - non-Hermitian eigenvalue problem
891: .  -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem 
892: -  -eps_pos_gen_non_hermitian - generalized non-Hermitian eigenvalue problem 
893:    with positive semi-definite B
894:     
895:    Notes:  
896:    Allowed values for the problem type are: Hermitian (EPS_HEP), non-Hermitian
897:    (EPS_NHEP), generalized Hermitian (EPS_GHEP), generalized non-Hermitian 
898:    (EPS_GNHEP), and generalized non-Hermitian with positive semi-definite B
899:    (EPS_PGNHEP).

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

908:    Level: beginner

910: .seealso: EPSSetOperators(), EPSSetType(), EPSGetProblemType(), EPSProblemType
911: @*/
912: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)
913: {
917:   switch (type) {
918:     case EPS_HEP:
919:       eps->isgeneralized = PETSC_FALSE;
920:       eps->ishermitian = PETSC_TRUE;
921:       eps->ispositive = PETSC_FALSE;
922:       break;
923:     case EPS_NHEP:
924:       eps->isgeneralized = PETSC_FALSE;
925:       eps->ishermitian = PETSC_FALSE;
926:       eps->ispositive = PETSC_FALSE;
927:       break;
928:     case EPS_GHEP:
929:       eps->isgeneralized = PETSC_TRUE;
930:       eps->ishermitian = PETSC_TRUE;
931:       eps->ispositive = PETSC_TRUE;
932:       break;
933:     case EPS_GNHEP:
934:       eps->isgeneralized = PETSC_TRUE;
935:       eps->ishermitian = PETSC_FALSE;
936:       eps->ispositive = PETSC_FALSE;
937:       break;
938:     case EPS_PGNHEP:
939:       eps->isgeneralized = PETSC_TRUE;
940:       eps->ishermitian = PETSC_FALSE;
941:       eps->ispositive = PETSC_TRUE;
942:       break;
943:     default:
944:       SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
945:   }
946:   eps->problem_type = type;
947:   return(0);
948: }

952: /*@C
953:    EPSGetProblemType - Gets the problem type from the EPS object.

955:    Not Collective

957:    Input Parameter:
958: .  eps - the eigensolver context 

960:    Output Parameter:
961: .  type - name of EPS problem type 

963:    Level: intermediate

965: .seealso: EPSSetProblemType(), EPSProblemType
966: @*/
967: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)
968: {
972:   *type = eps->problem_type;
973:   return(0);
974: }

978: /*@
979:    EPSSetExtraction - Specifies the type of extraction technique to be employed 
980:    by the eigensolver.

982:    Logically Collective on EPS

984:    Input Parameters:
985: +  eps  - the eigensolver context
986: -  extr - a known type of extraction

988:    Options Database Keys:
989: +  -eps_ritz - Rayleigh-Ritz extraction
990: .  -eps_harmonic - harmonic Ritz extraction
991: .  -eps_harmonic_relative - harmonic Ritz extraction relative to the eigenvalue
992: .  -eps_harmonic_right - harmonic Ritz extraction for rightmost eigenvalues
993: .  -eps_harmonic_largest - harmonic Ritz extraction for largest magnitude 
994:    (without target)
995: .  -eps_refined - refined Ritz extraction
996: -  -eps_refined_harmonic - refined harmonic Ritz extraction
997:     
998:    Notes:  
999:    Not all eigensolvers support all types of extraction. See the SLEPc
1000:    Users Manual for details.

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

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

1007:    Level: beginner

1009: .seealso: EPSSetTarget(), EPSGetExtraction(), EPSExtraction
1010: @*/
1011: PetscErrorCode EPSSetExtraction(EPS eps,EPSExtraction extr)
1012: {
1016:   eps->extraction = extr;
1017:   return(0);
1018: }

1022: /*@C
1023:    EPSGetExtraction - Gets the extraction type used by the EPS object.

1025:    Not Collective

1027:    Input Parameter:
1028: .  eps - the eigensolver context 

1030:    Output Parameter:
1031: .  extr - name of extraction type 

1033:    Level: intermediate

1035: .seealso: EPSSetExtraction(), EPSExtraction
1036: @*/
1037: PetscErrorCode EPSGetExtraction(EPS eps,EPSExtraction *extr)
1038: {
1042:   *extr = eps->extraction;
1043:   return(0);
1044: }

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

1052:    Logically Collective on EPS

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

1061:    Options Database Keys:
1062: +  -eps_balance <method> - the balancing method, where <method> is one of
1063:                            'none', 'oneside', 'twoside', or 'user'
1064: .  -eps_balance_its <its> - number of iterations
1065: -  -eps_balance_cutoff <cutoff> - cutoff value
1066:     
1067:    Notes:
1068:    When balancing is enabled, the solver works implicitly with matrix DAD^-1,
1069:    where D is an appropriate diagonal matrix. This improves the accuracy of
1070:    the computed results in some cases. See the SLEPc Users Manual for details.

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

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

1079:    The parameter 'its' is the number of iterations performed by the method. The
1080:    cutoff value is used only in the two-side variant. Use PETSC_IGNORE for an
1081:    argument that need not be changed. Use PETSC_DECIDE to assign a reasonably
1082:    good value.

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

1087:    Level: intermediate

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

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

1128:    Not Collective

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

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

1138:    Level: intermediate

1140:    Note:
1141:    The user can specify PETSC_NULL for any parameter that is not needed.

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

1157: /*@
1158:    EPSSetTrueResidual - Specifies if the solver must compute the true residual
1159:    explicitly or not.

1161:    Logically Collective on EPS

1163:    Input Parameters:
1164: +  eps     - the eigensolver context
1165: -  trueres - whether true residuals are required or not

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

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

1180:    Level: intermediate

1182: .seealso: EPSGetTrueResidual()
1183: @*/
1184: PetscErrorCode EPSSetTrueResidual(EPS eps,PetscBool trueres)
1185: {
1189:   eps->trueres = trueres;
1190:   return(0);
1191: }

1195: /*@C
1196:    EPSGetTrueResidual - Returns the flag indicating whether true
1197:    residuals must be computed explicitly or not.

1199:    Not Collective

1201:    Input Parameter:
1202: .  eps - the eigensolver context

1204:    Output Parameter:
1205: .  trueres - the returned flag

1207:    Level: intermediate

1209: .seealso: EPSSetTrueResidual()
1210: @*/
1211: PetscErrorCode EPSGetTrueResidual(EPS eps,PetscBool *trueres)
1212: {
1216:   *trueres = eps->trueres;
1217:   return(0);
1218: }

1222: /*@
1223:    EPSSetTrackAll - Specifies if the solver must compute the residual norm of all
1224:    approximate eigenpairs or not.

1226:    Logically Collective on EPS

1228:    Input Parameters:
1229: +  eps      - the eigensolver context
1230: -  trackall - whether to compute all residuals or not

1232:    Notes:
1233:    If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
1234:    the residual norm for each eigenpair approximation. Computing the residual is
1235:    usually an expensive operation and solvers commonly compute only the residual 
1236:    associated to the first unconverged eigenpair.

1238:    The options '-eps_monitor_all' and '-eps_monitor_draw_all' automatically
1239:    activate this option.

1241:    Level: intermediate

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

1256: /*@
1257:    EPSGetTrackAll - Returns the flag indicating whether all residual norms must
1258:    be computed or not.

1260:    Not Collective

1262:    Input Parameter:
1263: .  eps - the eigensolver context

1265:    Output Parameter:
1266: .  trackall - the returned flag

1268:    Level: intermediate

1270: .seealso: EPSSetTrackAll()
1271: @*/
1272: PetscErrorCode EPSGetTrackAll(EPS eps,PetscBool *trackall)
1273: {
1277:   *trackall = eps->trackall;
1278:   return(0);
1279: }

1283: /*@C
1284:    EPSSetOptionsPrefix - Sets the prefix used for searching for all 
1285:    EPS options in the database.

1287:    Logically Collective on EPS

1289:    Input Parameters:
1290: +  eps - the eigensolver context
1291: -  prefix - the prefix string to prepend to all EPS option requests

1293:    Notes:
1294:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1295:    The first character of all runtime options is AUTOMATICALLY the
1296:    hyphen.

1298:    For example, to distinguish between the runtime options for two
1299:    different EPS contexts, one could call
1300: .vb
1301:       EPSSetOptionsPrefix(eps1,"eig1_")
1302:       EPSSetOptionsPrefix(eps2,"eig2_")
1303: .ve

1305:    Level: advanced

1307: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
1308: @*/
1309: PetscErrorCode EPSSetOptionsPrefix(EPS eps,const char *prefix)
1310: {

1315:   if (!eps->OP) { EPSGetST(eps,&eps->OP); }
1316:   STSetOptionsPrefix(eps->OP,prefix);
1317:   PetscObjectSetOptionsPrefix((PetscObject)eps,prefix);
1318:   return(0);
1319: }
1320: 
1323: /*@C
1324:    EPSAppendOptionsPrefix - Appends to the prefix used for searching for all 
1325:    EPS options in the database.

1327:    Logically Collective on EPS

1329:    Input Parameters:
1330: +  eps - the eigensolver context
1331: -  prefix - the prefix string to prepend to all EPS option requests

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

1337:    Level: advanced

1339: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
1340: @*/
1341: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,const char *prefix)
1342: {

1347:   if (!eps->OP) { EPSGetST(eps,&eps->OP); }
1348:   STAppendOptionsPrefix(eps->OP,prefix);
1349:   PetscObjectAppendOptionsPrefix((PetscObject)eps,prefix);
1350:   return(0);
1351: }

1355: /*@C
1356:    EPSGetOptionsPrefix - Gets the prefix used for searching for all 
1357:    EPS options in the database.

1359:    Not Collective

1361:    Input Parameters:
1362: .  eps - the eigensolver context

1364:    Output Parameters:
1365: .  prefix - pointer to the prefix string used is returned

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

1370:    Level: advanced

1372: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
1373: @*/
1374: PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])
1375: {

1381:   PetscObjectGetOptionsPrefix((PetscObject)eps,prefix);
1382:   return(0);
1383: }