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-2012, 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 <slepc-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:     PetscOptionsBoolGroup("-eps_pos_gen_non_hermitian","generalized non-hermitian eigenvalue problem with positive semi-definite B","EPSSetProblemType",&flg);
 80:     if (flg) {EPSSetProblemType(eps,EPS_PGNHEP);}
 81:     PetscOptionsBoolGroupEnd("-eps_gen_indefinite","generalized hermitian-indefinite eigenvalue problem","EPSSetProblemType",&flg);
 82:     if (flg) {EPSSetProblemType(eps,EPS_GHIEP);}

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

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

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

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

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

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

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

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

227:   if (!eps->ip) { EPSGetIP(eps,&eps->ip); }
228:   IPSetFromOptions(eps->ip);
229:   if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
230:   DSSetFromOptions(eps->ds);
231:   STSetFromOptions(eps->OP);
232:   PetscRandomSetFromOptions(eps->rand);
233:   return(0);
234: }

238: /*@
239:    EPSGetTolerances - Gets the tolerance and maximum iteration count used
240:    by the EPS convergence tests. 

242:    Not Collective

244:    Input Parameter:
245: .  eps - the eigensolver context
246:   
247:    Output Parameters:
248: +  tol - the convergence tolerance
249: -  maxits - maximum number of iterations

251:    Notes:
252:    The user can specify PETSC_NULL for any parameter that is not needed.

254:    Level: intermediate

256: .seealso: EPSSetTolerances()
257: @*/
258: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,PetscInt *maxits)
259: {
262:   if (tol)    *tol    = eps->tol;
263:   if (maxits) *maxits = eps->max_it;
264:   return(0);
265: }

269: /*@
270:    EPSSetTolerances - Sets the tolerance and maximum iteration count used
271:    by the EPS convergence tests. 

273:    Logically Collective on EPS

275:    Input Parameters:
276: +  eps - the eigensolver context
277: .  tol - the convergence tolerance
278: -  maxits - maximum number of iterations to use

280:    Options Database Keys:
281: +  -eps_tol <tol> - Sets the convergence tolerance
282: -  -eps_max_it <maxits> - Sets the maximum number of iterations allowed

284:    Notes:
285:    Use PETSC_IGNORE for an argument that need not be changed.

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

290:    Level: intermediate

292: .seealso: EPSGetTolerances()
293: @*/
294: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,PetscInt maxits)
295: {
300:   if (tol != PETSC_IGNORE) {
301:     if (tol == PETSC_DEFAULT) {
302:       eps->tol = PETSC_DEFAULT;
303:     } else {
304:       if (tol < 0.0) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
305:       eps->tol = tol;
306:     }
307:   }
308:   if (maxits != PETSC_IGNORE) {
309:     if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
310:       eps->max_it = 0;
311:       eps->setupcalled = 0;
312:     } else {
313:       if (maxits < 0) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
314:       eps->max_it = maxits;
315:     }
316:   }
317:   return(0);
318: }

322: /*@
323:    EPSGetDimensions - Gets the number of eigenvalues to compute
324:    and the dimension of the subspace.

326:    Not Collective

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

336:    Notes:
337:    The user can specify PETSC_NULL for any parameter that is not needed.

339:    Level: intermediate

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

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

359:    Logically Collective on EPS

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

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

372:    Notes:
373:    Use PETSC_IGNORE to retain the previous value of any parameter.

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

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

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

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

394:    Level: intermediate

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

432: /*@
433:    EPSSetWhichEigenpairs - Specifies which portion of the spectrum is 
434:    to be sought.

436:    Logically Collective on EPS

438:    Input Parameters:
439: +  eps   - eigensolver context obtained from EPSCreate()
440: -  which - the portion of the spectrum to be sought

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

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

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

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

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

486:    Level: intermediate

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

527: /*@C
528:    EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be 
529:    sought.

531:    Not Collective

533:    Input Parameter:
534: .  eps - eigensolver context obtained from EPSCreate()

536:    Output Parameter:
537: .  which - the portion of the spectrum to be sought

539:    Notes:
540:    See EPSSetWhichEigenpairs() for possible values of 'which'.

542:    Level: intermediate

544: .seealso: EPSSetWhichEigenpairs(), EPSWhich
545: @*/
546: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)
547: {
551:   *which = eps->which;
552:   return(0);
553: }

557: /*@
558:    EPSSetLeftVectorsWanted - Specifies which eigenvectors are required.

560:    Logically Collective on EPS

562:    Input Parameters:
563: +  eps      - the eigensolver context
564: -  leftvecs - whether left eigenvectors are required or not

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

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

574:    Level: intermediate

576: .seealso: EPSGetLeftVectorsWanted(), EPSGetEigenvectorLeft()
577: @*/
578: PetscErrorCode EPSSetLeftVectorsWanted(EPS eps,PetscBool leftvecs)
579: {
583:   if (eps->leftvecs != leftvecs) {
584:     eps->leftvecs = leftvecs;
585:     eps->setupcalled = 0;
586:   }
587:   return(0);
588: }

592: /*@
593:    EPSGetLeftVectorsWanted - Returns the flag indicating whether left 
594:    eigenvectors are required or not.

596:    Not Collective

598:    Input Parameter:
599: .  eps - the eigensolver context

601:    Output Parameter:
602: .  leftvecs - the returned flag

604:    Level: intermediate

606: .seealso: EPSSetLeftVectorsWanted(), EPSWhich
607: @*/
608: PetscErrorCode EPSGetLeftVectorsWanted(EPS eps,PetscBool *leftvecs)
609: {
613:   *leftvecs = eps->leftvecs;
614:   return(0);
615: }

619: /*@
620:    EPSSetMatrixNorms - Gives the reference values of the matrix norms
621:    and specifies whether these values should be improved adaptively.

623:    Logically Collective on EPS

625:    Input Parameters:
626: +  eps      - the eigensolver context
627: .  nrma     - a reference value for the norm of matrix A
628: .  nrmb     - a reference value for the norm of matrix B
629: -  adaptive - whether matrix norms are improved adaptively

631:    Options Database Keys:
632: +  -eps_norm_a <nrma> - norm of A
633: .  -eps_norm_b <nrma> - norm of B
634: -  -eps_norms_adaptive <boolean> - Sets/resets the boolean flag 'adaptive'

636:    Notes:
637:    If the user sets adaptive=PETSC_FALSE then the solver uses the values
638:    of nrma and nrmb for the matrix norms, and these values do not change
639:    throughout the iteration.

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

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

649:    Level: intermediate

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

691: /*@
692:    EPSGetMatrixNorms - Returns the value of the matrix norms (either set
693:    by the user or estimated by the solver) and the flag indicating whether
694:    the norms are being adaptively improved.

696:    Not Collective

698:    Input Parameter:
699: .  eps - the eigensolver context

701:    Output Parameters:
702: +  nrma     - the norm of matrix A
703: .  nrmb     - the norm of matrix B
704: -  adaptive - whether matrix norms are improved adaptively

706:    Level: intermediate

708: .seealso: EPSSetMatrixNorms()
709: @*/
710: PetscErrorCode EPSGetMatrixNorms(EPS eps,PetscReal *nrma,PetscReal *nrmb,PetscBool *adaptive)
711: {
714:   if (nrma) *nrma = eps->nrma;
715:   if (nrmb) *nrmb = eps->nrmb;
716:   if (adaptive) *adaptive = eps->adaptive;
717:   return(0);
718: }

722: /*@C
723:    EPSSetEigenvalueComparison - Specifies the eigenvalue comparison function
724:    when EPSSetWhichEigenpairs() is set to EPS_WHICH_USER.

726:    Logically Collective on EPS

728:    Input Parameters:
729: +  eps  - eigensolver context obtained from EPSCreate()
730: .  func - a pointer to the comparison function
731: -  ctx  - a context pointer (the last parameter to the comparison function)

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

736: +   ar     - real part of the 1st eigenvalue
737: .   ai     - imaginary part of the 1st eigenvalue
738: .   br     - real part of the 2nd eigenvalue
739: .   bi     - imaginary part of the 2nd eigenvalue
740: .   res    - result of comparison
741: -   ctx    - optional context, as set by EPSSetEigenvalueComparison()

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

749:    Level: advanced

751: .seealso: EPSSetWhichEigenpairs(), EPSSortEigenvalues(), EPSWhich
752: @*/
753: PetscErrorCode EPSSetEigenvalueComparison(EPS eps,PetscErrorCode (*func)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)
754: {
757:   eps->which_func = func;
758:   eps->which_ctx  = ctx;
759:   eps->which      = EPS_WHICH_USER;
760:   return(0);
761: }

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

770:    Logically Collective on EPS

772:    Input Parameters:
773: +  eps  - eigensolver context obtained from EPSCreate()
774: .  func - a pointer to the evaluation function
775: -  ctx  - a context pointer (the last parameter to the evaluation function)

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

780: +   er     - real part of the current eigenvalue approximation
781: .   ei     - imaginary part of the current eigenvalue approximation
782: .   xr     - real part of the current eigenvector approximation
783: .   xi     - imaginary part of the current eigenvector approximation
784: .   rr     - result of evaluation (real part)
785: .   ri     - result of evaluation (imaginary part)
786: -   ctx    - optional context, as set by EPSSetArbitrarySelection()

788:    Note:
789:    This provides a mechanism to select eigenpairs by evaluating a user-defined
790:    function. When a function has been provided, the default selection based on
791:    sorting the eigenvalues is replaced by the sorting of the results of this
792:    function (with the same sorting criterion given in EPSSortEigenvalues()).

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

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

804:    Level: advanced

806: .seealso: EPSSetWhichEigenpairs(), EPSSortEigenvalues()
807: @*/
808: PetscErrorCode EPSSetArbitrarySelection(EPS eps,PetscErrorCode (*func)(PetscScalar,PetscScalar,Vec,Vec,PetscScalar*,PetscScalar*,void*),void* ctx)
809: {
812:   eps->arbit_func = func;
813:   eps->arbit_ctx  = ctx;
814:   eps->setupcalled = 0;
815:   return(0);
816: }

820: /*@C
821:    EPSSetConvergenceTestFunction - Sets a function to compute the error estimate
822:    used in the convergence test.

824:    Logically Collective on EPS

826:    Input Parameters:
827: +  eps  - eigensolver context obtained from EPSCreate()
828: .  func - a pointer to the convergence test function
829: -  ctx  - a context pointer (the last parameter to the convergence test function)

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

834: +   eps    - eigensolver context obtained from EPSCreate()
835: .   eigr   - real part of the eigenvalue
836: .   eigi   - imaginary part of the eigenvalue
837: .   res    - residual norm associated to the eigenpair
838: .   errest - (output) computed error estimate
839: -   ctx    - optional context, as set by EPSSetConvergenceTest()

841:    Note:
842:    If the error estimate returned by the convergence test function is less than
843:    the tolerance, then the eigenvalue is accepted as converged.

845:    Level: advanced

847: .seealso: EPSSetConvergenceTest(),EPSSetTolerances()
848: @*/
849: extern PetscErrorCode EPSSetConvergenceTestFunction(EPS eps,PetscErrorCode (*func)(EPS,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx)
850: {
853:   eps->conv_func = func;
854:   eps->conv_ctx = ctx;
855:   if (func == EPSConvergedEigRelative) eps->conv = EPS_CONV_EIG;
856:   else if (func == EPSConvergedNormRelative) eps->conv = EPS_CONV_NORM;
857:   else if (func == EPSConvergedAbsolute) eps->conv = EPS_CONV_ABS;
858:   else eps->conv = EPS_CONV_USER;
859:   return(0);
860: }

864: /*@
865:    EPSSetConvergenceTest - Specifies how to compute the error estimate
866:    used in the convergence test.

868:    Logically Collective on EPS

870:    Input Parameters:
871: +  eps   - eigensolver context obtained from EPSCreate()
872: -  conv  - the type of convergence test

874:    Options Database Keys:
875: +  -eps_conv_abs - Sets the absolute convergence test
876: .  -eps_conv_eig - Sets the convergence test relative to the eigenvalue
877: -  -eps_conv_norm - Sets the convergence test relative to the matrix norms

879:    Note:
880:    The parameter 'conv' can have one of these values
881: +     EPS_CONV_ABS - absolute error ||r||
882: .     EPS_CONV_EIG - error relative to the eigenvalue l, ||r||/|l|
883: .     EPS_CONV_NORM - error relative to the matrix norms, ||r||/(||A||+|l|*||B||)
884: -     EPS_CONV_USER - function set by EPSSetConvergenceTestFunction() 

886:    Level: intermediate

888: .seealso: EPSGetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSConv
889: @*/
890: PetscErrorCode EPSSetConvergenceTest(EPS eps,EPSConv conv)
891: {
895:   switch(conv) {
896:     case EPS_CONV_EIG: eps->conv_func = EPSConvergedEigRelative; break;
897:     case EPS_CONV_NORM: eps->conv_func = EPSConvergedNormRelative; break;
898:     case EPS_CONV_ABS: eps->conv_func = EPSConvergedAbsolute; break;
899:     case EPS_CONV_USER: break;
900:     default:
901:       SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
902:   }
903:   eps->conv = conv;
904:   return(0);
905: }

909: /*@
910:    EPSGetConvergenceTest - Gets the method used to compute the error estimate
911:    used in the convergence test.

913:    Not Collective

915:    Input Parameters:
916: .  eps   - eigensolver context obtained from EPSCreate()

918:    Output Parameters:
919: .  conv  - the type of convergence test

921:    Level: intermediate

923: .seealso: EPSSetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSConv
924: @*/
925: PetscErrorCode EPSGetConvergenceTest(EPS eps,EPSConv *conv)
926: {
930:   *conv = eps->conv;
931:   return(0);
932: }

936: /*@
937:    EPSSetProblemType - Specifies the type of the eigenvalue problem.

939:    Logically Collective on EPS

941:    Input Parameters:
942: +  eps      - the eigensolver context
943: -  type     - a known type of eigenvalue problem 

945:    Options Database Keys:
946: +  -eps_hermitian - Hermitian eigenvalue problem
947: .  -eps_gen_hermitian - generalized Hermitian eigenvalue problem
948: .  -eps_non_hermitian - non-Hermitian eigenvalue problem
949: .  -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem 
950: -  -eps_pos_gen_non_hermitian - generalized non-Hermitian eigenvalue problem 
951:    with positive semi-definite B
952:     
953:    Notes:  
954:    Allowed values for the problem type are: Hermitian (EPS_HEP), non-Hermitian
955:    (EPS_NHEP), generalized Hermitian (EPS_GHEP), generalized non-Hermitian 
956:    (EPS_GNHEP), generalized non-Hermitian with positive semi-definite B
957:    (EPS_PGNHEP), and generalized Hermitian-indefinite (EPS_GHIEP).

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

966:    Level: beginner

968: .seealso: EPSSetOperators(), EPSSetType(), EPSGetProblemType(), EPSProblemType
969: @*/
970: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)
971: {
975:   switch (type) {
976:     case EPS_HEP:
977:       eps->isgeneralized = PETSC_FALSE;
978:       eps->ishermitian = PETSC_TRUE;
979:       eps->ispositive = PETSC_FALSE;
980:       break;
981:     case EPS_NHEP:
982:       eps->isgeneralized = PETSC_FALSE;
983:       eps->ishermitian = PETSC_FALSE;
984:       eps->ispositive = PETSC_FALSE;
985:       break;
986:     case EPS_GHEP:
987:       eps->isgeneralized = PETSC_TRUE;
988:       eps->ishermitian = PETSC_TRUE;
989:       eps->ispositive = PETSC_TRUE;
990:       break;
991:     case EPS_GNHEP:
992:       eps->isgeneralized = PETSC_TRUE;
993:       eps->ishermitian = PETSC_FALSE;
994:       eps->ispositive = PETSC_FALSE;
995:       break;
996:     case EPS_PGNHEP:
997:       eps->isgeneralized = PETSC_TRUE;
998:       eps->ishermitian = PETSC_FALSE;
999:       eps->ispositive = PETSC_TRUE;
1000:       break;
1001:     case EPS_GHIEP:
1002:       eps->isgeneralized = PETSC_TRUE;
1003:       eps->ishermitian = PETSC_TRUE;
1004:       eps->ispositive = PETSC_FALSE;
1005:       break;
1006:     default:
1007:       SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
1008:   }
1009:   eps->problem_type = type;
1010:   return(0);
1011: }

1015: /*@C
1016:    EPSGetProblemType - Gets the problem type from the EPS object.

1018:    Not Collective

1020:    Input Parameter:
1021: .  eps - the eigensolver context 

1023:    Output Parameter:
1024: .  type - name of EPS problem type 

1026:    Level: intermediate

1028: .seealso: EPSSetProblemType(), EPSProblemType
1029: @*/
1030: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)
1031: {
1035:   *type = eps->problem_type;
1036:   return(0);
1037: }

1041: /*@
1042:    EPSSetExtraction - Specifies the type of extraction technique to be employed 
1043:    by the eigensolver.

1045:    Logically Collective on EPS

1047:    Input Parameters:
1048: +  eps  - the eigensolver context
1049: -  extr - a known type of extraction

1051:    Options Database Keys:
1052: +  -eps_ritz - Rayleigh-Ritz extraction
1053: .  -eps_harmonic - harmonic Ritz extraction
1054: .  -eps_harmonic_relative - harmonic Ritz extraction relative to the eigenvalue
1055: .  -eps_harmonic_right - harmonic Ritz extraction for rightmost eigenvalues
1056: .  -eps_harmonic_largest - harmonic Ritz extraction for largest magnitude 
1057:    (without target)
1058: .  -eps_refined - refined Ritz extraction
1059: -  -eps_refined_harmonic - refined harmonic Ritz extraction
1060:     
1061:    Notes:  
1062:    Not all eigensolvers support all types of extraction. See the SLEPc
1063:    Users Manual for details.

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

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

1070:    Level: beginner

1072: .seealso: EPSSetTarget(), EPSGetExtraction(), EPSExtraction
1073: @*/
1074: PetscErrorCode EPSSetExtraction(EPS eps,EPSExtraction extr)
1075: {
1079:   eps->extraction = extr;
1080:   return(0);
1081: }

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

1088:    Not Collective

1090:    Input Parameter:
1091: .  eps - the eigensolver context 

1093:    Output Parameter:
1094: .  extr - name of extraction type 

1096:    Level: intermediate

1098: .seealso: EPSSetExtraction(), EPSExtraction
1099: @*/
1100: PetscErrorCode EPSGetExtraction(EPS eps,EPSExtraction *extr)
1101: {
1105:   *extr = eps->extraction;
1106:   return(0);
1107: }

1111: /*@
1112:    EPSSetBalance - Specifies the balancing technique to be employed by the
1113:    eigensolver, and some parameters associated to it.

1115:    Logically Collective on EPS

1117:    Input Parameters:
1118: +  eps    - the eigensolver context
1119: .  bal    - the balancing method, one of EPS_BALANCE_NONE, EPS_BALANCE_ONESIDE,
1120:             EPS_BALANCE_TWOSIDE, or EPS_BALANCE_USER
1121: .  its    - number of iterations of the balancing algorithm
1122: -  cutoff - cutoff value

1124:    Options Database Keys:
1125: +  -eps_balance <method> - the balancing method, where <method> is one of
1126:                            'none', 'oneside', 'twoside', or 'user'
1127: .  -eps_balance_its <its> - number of iterations
1128: -  -eps_balance_cutoff <cutoff> - cutoff value
1129:     
1130:    Notes:
1131:    When balancing is enabled, the solver works implicitly with matrix DAD^-1,
1132:    where D is an appropriate diagonal matrix. This improves the accuracy of
1133:    the computed results in some cases. See the SLEPc Users Manual for details.

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

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

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

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

1150:    Level: intermediate

1152: .seealso: EPSGetBalance(), EPSBalance, STSetBalanceMatrix()
1153: @*/
1154: PetscErrorCode EPSSetBalance(EPS eps,EPSBalance bal,PetscInt its,PetscReal cutoff)
1155: {
1161:   if (bal!=PETSC_IGNORE) {
1162:     if (bal==PETSC_DECIDE || bal==PETSC_DEFAULT) eps->balance = EPS_BALANCE_TWOSIDE;
1163:     else switch (bal) {
1164:       case EPS_BALANCE_NONE:
1165:       case EPS_BALANCE_ONESIDE:
1166:       case EPS_BALANCE_TWOSIDE:
1167:       case EPS_BALANCE_USER:
1168:         eps->balance = bal;
1169:         break;
1170:       default:
1171:         SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid value of argument 'bal'");
1172:     }
1173:   }
1174:   if (its!=PETSC_IGNORE) {
1175:     if (its==PETSC_DECIDE || its==PETSC_DEFAULT) eps->balance_its = 5;
1176:     else eps->balance_its = its;
1177:   }
1178:   if (cutoff!=PETSC_IGNORE) {
1179:     if (cutoff==PETSC_DECIDE || cutoff==PETSC_DEFAULT) eps->balance_cutoff = 1e-8;
1180:     else eps->balance_cutoff = cutoff;
1181:   }
1182:   return(0);
1183: }

1187: /*@
1188:    EPSGetBalance - Gets the balancing type used by the EPS object, and the associated
1189:    parameters.

1191:    Not Collective

1193:    Input Parameter:
1194: .  eps - the eigensolver context 

1196:    Output Parameters:
1197: +  bal    - the balancing method
1198: .  its    - number of iterations of the balancing algorithm
1199: -  cutoff - cutoff value

1201:    Level: intermediate

1203:    Note:
1204:    The user can specify PETSC_NULL for any parameter that is not needed.

1206: .seealso: EPSSetBalance(), EPSBalance
1207: @*/
1208: PetscErrorCode EPSGetBalance(EPS eps,EPSBalance *bal,PetscInt *its,PetscReal *cutoff)
1209: {
1212:   if (bal)    *bal = eps->balance;
1213:   if (its)    *its = eps->balance_its;
1214:   if (cutoff) *cutoff = eps->balance_cutoff;
1215:   return(0);
1216: }

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

1224:    Logically Collective on EPS

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

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

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

1243:    Level: intermediate

1245: .seealso: EPSGetTrueResidual()
1246: @*/
1247: PetscErrorCode EPSSetTrueResidual(EPS eps,PetscBool trueres)
1248: {
1252:   eps->trueres = trueres;
1253:   return(0);
1254: }

1258: /*@C
1259:    EPSGetTrueResidual - Returns the flag indicating whether true
1260:    residuals must be computed explicitly or not.

1262:    Not Collective

1264:    Input Parameter:
1265: .  eps - the eigensolver context

1267:    Output Parameter:
1268: .  trueres - the returned flag

1270:    Level: intermediate

1272: .seealso: EPSSetTrueResidual()
1273: @*/
1274: PetscErrorCode EPSGetTrueResidual(EPS eps,PetscBool *trueres)
1275: {
1279:   *trueres = eps->trueres;
1280:   return(0);
1281: }

1285: /*@
1286:    EPSSetTrackAll - Specifies if the solver must compute the residual norm of all
1287:    approximate eigenpairs or not.

1289:    Logically Collective on EPS

1291:    Input Parameters:
1292: +  eps      - the eigensolver context
1293: -  trackall - whether to compute all residuals or not

1295:    Notes:
1296:    If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
1297:    the residual norm for each eigenpair approximation. Computing the residual is
1298:    usually an expensive operation and solvers commonly compute only the residual 
1299:    associated to the first unconverged eigenpair.

1301:    The options '-eps_monitor_all' and '-eps_monitor_draw_all' automatically
1302:    activate this option.

1304:    Level: intermediate

1306: .seealso: EPSGetTrackAll()
1307: @*/
1308: PetscErrorCode EPSSetTrackAll(EPS eps,PetscBool trackall)
1309: {
1313:   eps->trackall = trackall;
1314:   return(0);
1315: }

1319: /*@
1320:    EPSGetTrackAll - Returns the flag indicating whether all residual norms must
1321:    be computed or not.

1323:    Not Collective

1325:    Input Parameter:
1326: .  eps - the eigensolver context

1328:    Output Parameter:
1329: .  trackall - the returned flag

1331:    Level: intermediate

1333: .seealso: EPSSetTrackAll()
1334: @*/
1335: PetscErrorCode EPSGetTrackAll(EPS eps,PetscBool *trackall)
1336: {
1340:   *trackall = eps->trackall;
1341:   return(0);
1342: }

1346: /*@C
1347:    EPSSetOptionsPrefix - Sets the prefix used for searching for all 
1348:    EPS options in the database.

1350:    Logically Collective on EPS

1352:    Input Parameters:
1353: +  eps - the eigensolver context
1354: -  prefix - the prefix string to prepend to all EPS option requests

1356:    Notes:
1357:    A hyphen (-) must NOT be given at the beginning of the prefix name.
1358:    The first character of all runtime options is AUTOMATICALLY the
1359:    hyphen.

1361:    For example, to distinguish between the runtime options for two
1362:    different EPS contexts, one could call
1363: .vb
1364:       EPSSetOptionsPrefix(eps1,"eig1_")
1365:       EPSSetOptionsPrefix(eps2,"eig2_")
1366: .ve

1368:    Level: advanced

1370: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
1371: @*/
1372: PetscErrorCode EPSSetOptionsPrefix(EPS eps,const char *prefix)
1373: {

1378:   if (!eps->OP) { EPSGetST(eps,&eps->OP); }
1379:   STSetOptionsPrefix(eps->OP,prefix);
1380:   if (!eps->ip) { EPSGetIP(eps,&eps->ip); }
1381:   IPSetOptionsPrefix(eps->ip,prefix);
1382:   if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
1383:   DSSetOptionsPrefix(eps->ds,prefix);
1384:   PetscObjectSetOptionsPrefix((PetscObject)eps,prefix);
1385:   return(0);
1386: }
1387: 
1390: /*@C
1391:    EPSAppendOptionsPrefix - Appends to the prefix used for searching for all 
1392:    EPS options in the database.

1394:    Logically Collective on EPS

1396:    Input Parameters:
1397: +  eps - the eigensolver context
1398: -  prefix - the prefix string to prepend to all EPS option requests

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

1404:    Level: advanced

1406: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
1407: @*/
1408: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,const char *prefix)
1409: {

1414:   if (!eps->OP) { EPSGetST(eps,&eps->OP); }
1415:   STAppendOptionsPrefix(eps->OP,prefix);
1416:   if (!eps->ip) { EPSGetIP(eps,&eps->ip); }
1417:   IPSetOptionsPrefix(eps->ip,prefix);
1418:   if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
1419:   DSSetOptionsPrefix(eps->ds,prefix);
1420:   PetscObjectAppendOptionsPrefix((PetscObject)eps,prefix);
1421:   return(0);
1422: }

1426: /*@C
1427:    EPSGetOptionsPrefix - Gets the prefix used for searching for all 
1428:    EPS options in the database.

1430:    Not Collective

1432:    Input Parameters:
1433: .  eps - the eigensolver context

1435:    Output Parameters:
1436: .  prefix - pointer to the prefix string used is returned

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

1441:    Level: advanced

1443: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
1444: @*/
1445: PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])
1446: {

1452:   PetscObjectGetOptionsPrefix((PetscObject)eps,prefix);
1453:   return(0);
1454: }